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SUMMARY 


Presented in this report is a comprehensive review of the 
existing technical literature and a design summary of stall 
characteristics applicable to light straight wing aircraft. 

These characteristics are obtained with the aid of a digital 
computer program which utilizes the most up to date analytical 
methods employing lifting line theory and the available exper- 
imental test data for wing section characteristics. The computer 
results are presented in the form of stall charts suitable for 
preliminary design purposes. Based on the extensive parametric 
study covering a total of 331 different aircraft configurations, 
it can be concluded that in modern airplane design satisfactory 
stalling characteristics can be readily built in with no apprecia- 
ble loss in airplane performance or handling qualities. A proper 
combination of wing taper, twist and type of airfoil sections 
with minor post-design fixes, if required, can in most cases 
provide satisfactory wing stall characteristics . 
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The work was performed by the Scientific Systems Division 
(SSD) of the Dynasciences Corporation, Blue Bell, Pennsylvania , 
for the National Aeronautics and Space Administration (NASA), 
Langley Research Center, Hampton, Virginia, under contract number 
NAS 1-8389 during the period from July 1968 through September 
1969. 


The NASA technical representatives were Mr. Robert T. Taylor 
and Mr. William J. Alford, Jr. The contributions of the NASA 
technical personnel to this work are gratefully acknowledged. 
Acknowledgement is also extended to NASA computer personnel, 
especially Mrs. Belinda Adams, for their support in this program. 

Messrs. James C. Sivells and Hartley A. Soule were special 
technical consultants on this project and Mr. Ron Anton was 
computer consultant . 


v 





CONTENTS 


Page 

SUMMARY iii 

FOREWORD v 

LIST OF ILLUSTRATIONS viii 

LIST OF TABLES xi 

LIST OF SYMBOLS xiii 

SECTION 1 INTRODUCTION 1 

SECTION 2 BASIC CONSIDERATIONS OF AIRPLANE 

STALLING 4 

SECTION 3 THEORETICAL ANALYSIS 21 

SECTION 4 COMPUTER PROGRAM 42 

SECTION 5 PARAMETRIC INVESTIGATION 75 

131 

SECTION 6 SCALE MODEL WIND TUNNEL TESTING 

SECTION 7 DESIGN PROCEDURES 136 

SECTION 8 CONCLUSIONS AND RECOMMENDATIONS 142 

SECTION 9 REFERENCES 144 

APPENDIX A INTERNAL LISTING OF THE COMPUTER 

PROGRAM 

APPENDIX B AVAILABILITY OF THE COMPUTER PROGRAM 189 


BIBLIOGRAPHY 


190 


ILLUSTRATIONS 


Figure Page 


1 Representative Lift Curve (Reproduced from 

Reference 13) 10 

2 The Low-Speed Stalling Characteristics of Airfoil 
Sections Correlated With. Reynolds Number and 

the Upper- Surf ace Ordinates of the Airfoil Sections 
at the 0. 0125-Chord Station 12 

3 Wing Leading Edge Mof if icat ions for Controlling 

Wing Stall 20 

4 Definition of Parameters for Transformation of 

Wing-Body Combination 26 


5 Typical Load Distributions for Obtaining Factors 

for Altering Two-Dimensional Data 


6 Illustration of Method for Correcting Two- 

Dimensional Section Data 39 

7 Extrapolations of Lift Curve Slopes at Low 

Reynolds Number 51 

8 Variation of Section Lift-Curve Slope with Thick- 
ness-Chord Ratio at Constant Reynolds Number 

NACA 644 Sections 54 

9 Corrected Lift Curves for NACA 64-421 Airfoil 

at Low Reynolds Numbers 55 


10 Method of Tabulation of Section Characteristics .. 57 

11 Schematic Representation of Section Data Storage 


in the Computer 59 

12 Nomenclature for Developing Interpolation Formulae 61 

13 Computer Program Block Diagram 62 

14 Schematic Representation of the Computer Input 

Cards 63 

15 Experimental and Calculated Characteristics for 

a Wing of Aspect Ratio 8.04 68 


viii 



Figure 


Experimental and Calculated Characteristics 

for a Wing of Aspect Ratio 10.05 70 

Experimental and Calculated Characteristics 

for a Wing of Aspect Ratio 12.06 72 

Experimental and Calculated Characteristics 

for Wing with 60% Flap; Aspect Ratio 9.02; Taper 

Ratio 0.4; Washout 2° 74 

Typical Lift Distributions Along Wing Span 81 

Variation of c l max with Reynolds Number and 
Thickness -Chord Ratio . 83 

Variation of 0L max with Reynolds Number and Taper 
Ratio • • • 


Effect of Aspect Ratio on Stall Margin 
Distribution 


Effect of Aspect Ratio on Wing Stall Pattern .... 91 

Effect of Aspect Ratio and Taper Ratio on C L m ax. . 94 

Increment of Induced-Drag Coefficient Due to 
Washout, 230 Series Airfoil Section 98 


Effect of Root Thickness-Chord Ratio on Stall 
Margin Distribution ^01- 

Effect of Root Thickness -Chord Ratio on Wing 

Stall Boundaries 104 

Effect of Root Thickness -Chord Ratio on 0L m a X . . . 107 

Effect of Tip Thickness-Chord Ratio on Stall 
Margin Distribution 1- 11- 

Effect of Tip Thickness-Chord Ratio on^L-max ... 114 

Effect of Reynolds Number on Stall Margin 
Distribution 117 


Effect of Reynolds Number on Wing Stall 
Pattern 


ix 



Figure Page 

33 Effect of Reynolds Number on C -fmay 123 

34 Effect of Wing Camber on Stall Margin 

Distribution 127 

35 Effect of Fuselage 128 

36 Effect of the Span of a 20% Chord Split Flap 

on the Wing Stalling Characteristics 130 

37 Variation of Wing Maximum Lift Coefficient 

with. Stall Speed 138 


x 



TABLES 


Tables Page 

I Airfoil Section Data Available for Use with 

the Computer Program 56 

II Typical Computer Output 65 

III Summary of Configurations Studied 76 




SYMBOLS 


non-dimensional fuselage semiheight 
fuselage semiheight , ft . 
wing aspect ratio, 

coefficients in trigonometric series 

non-dimensional average distance of ppint on 
wing from fuselage cross-section focii 

average distance of point on wing from fuselage 
cross-section focii, ft. 

section lift-curve slope, per degree 

non-dimensional fuselage semiwidth 

fuselage semiwidth, f t . 

wing span, ft. 

flap span, ft. 

total wing drag coefficient, 

wing profile drag coefficient,-^ 

qS 

wing induced drag coefficient,-^- 
wing lift coefficient, 

qb 

wing pitching moment coefficient, 
total section drag coefficient 
section induced drag coefficient 
section profile drag coefficient 
section lift coefficient, 


C it 
C/o 


^ mox 
( C * max ) 0 

C i* 

C 'l 


C ' 2 


C/ a 

c/8 



c 

c 1 


C R 

C T 

D 

Do 
D i 


E 

e 


design lift coefficient 

two-dimensional, uncorrected value of 
lift coefficient 

maximum section lift coefficient 

maximum two-dimensional section lift 
coefficient 

section lift coefficient at end of flap 

section lift coefficient for that part 
of lift distribution involving no 
discontinuity in angle of attack 

section lift coefficient for part of 
lift distribution due to discontinuity 
in angle of attack 

section lift-curve slope, per degree 

section lift coefficient with deflected 
flaps, calculated assuming linear lift 
curve s 

section pitching moment coefficient 

section pitching moment coefficient 
about quarter-chord point 

wing chord at any spanwise station, ft. 

wing mean aerodynamic chord, ft. 

wing root chord, ft. 

wing tip chord, ft. 

total wing drag, lbs. 

wing profile drag, lbs. 

wing induced drag, lbs. 

edge velocity factor 

non-dimensional eccentricity of fuselage 
cross-section 


xiv 


e 1 

F 

FF 

Gmk 

^mk 

H 

K 


K ij 


s 

T 

t 

tr 

X 

Y 

Y 0 

y,y 


eccentricity of fuselage cross-section, ft. 

factor used in altering two-dimensional 
lift curves 

C / max 

value of jz7 r at y = y* taken at the flap 

side ofy* (C/mQx) ° 

coefficient of transpose of matrix G'mk 
coefficient of a matrix 

non-dimensional wing height above fuselage 
centerline 

section camber level measured in terms of 
section design lift coefficient 

coefficient of the inverse of matrix ^mk 

wing lift, lbs. 

section lift, lbs. 

wing pitching moment, lbs. ft. 

2 

dynamic pressure, lbs/ft. 

Reynolds number 

Reynolds number based on mean aerodynamic chord 
gross wing area, sq. ft. 
thickness factor 

wing section maximum thickness, ft. 

maximum wing root thickness, ft. 

coordinate parallel to fuselage centerline, 
positive forward, ft. 

non-dimensional spanwise distance 

non-dimensional coordinate of wing-fuselage 
junction 

spanwise coordinate, ft. 


x v 



A DESIGN SUMMARY OF STALL CHARACTERISTICS 
OF STRAIGHT WING AIRCRAFT 

By M. A. McVeigh, E. Kisielowski 
DYNASC IENCE S CORPORATION 

SECTION I 
INTRODUCTION 

Federal and military aviation regulations require that all 
aircraft possess stalling characteristics which comply with the 
established specifications. The final evaluation of an airplane 
is made in flight tests in which all factors influencing its 
stalling characteristics are integrated. At this stage, however, 
it may be too late or very expensive to make any alterations 
necessary to meet the established requirements. 

It is highly desirable therefore, to predict the maximum lift 
and stalling characteristics of an aircraft at an early design 
stage and, on the basis of the predicted characteristics, to 
incorporate such modifications in the design as may give promise 
of satisfying the FAA or military stall requirements. With this 
end in mind, extensive research effort directed toward a compre- 
hensive understanding of the aerodynamics of unswept wings has 
been accomplished in prior years. This research effort was dir- 
ected along lines of theoretical work, wind tunnel experiments and 
flight tests. 

On the theoretical front substantial progress has been made. 
Prior to the World War II an adaptation of the lifting line theory 
in which the wing airfoil sections were assumed to possess a linear 
variation of lift with angle-of-attack was the only available 
method for predicting the maximum lift and stalling characteristics 
of unswept wings. In the region of maximum lift, however, the 
section lift curves are usually quite nonlinear. Consequently, 
predictions based on the method were quite unreliable in the max- 
imum lift range. Shortly after the close of World War II Sivells 
(References 1 and 2) presented a method for calculating unswept 
wing characteristics by lifting line theory utilizing nonlinear 
airfoil section lift data. The values of maximum lift and the 
point of initial stall predicted by this approach agree well with 
experimental observations, and the method is regarded as a very 
effective design tool. Other major theoretical contributions were 
made by Theodorsen (Reference 3) and by Multhopp (Reference 4). 
These references provide valuable methods for developing improved 
airfoil sections and for determining wing-fuselage interference 
effects. 

Wind tunnel experiments have provided extensive information 
in such important areas as: two-dimensional airfoil section 


characteristics (References 5 and 6); the aerodynamic charac- 
teristics of complete wings at high Reynolds number (Reference 
7); the characteristics of a large number of complete airplanes 
(Reference 8); and compressibility effects on maximum lift (Ref- 
erence 9). Certain portions of the wind tunnel test data have 
direct design applications, other portions provide a basis on 
which to evaluate the applicability of new theories, and still 
other portions provide information on which to base empirical 
design guidance. 

The flight test results provide a comprehensive evaluation 
of the integrated effect of all the factors that contribute to 
the airplane flying qualities. Measurement of the flight charac- 
teristics of many airplanes of different types has permitted a 
definition of those characteristics which provide for good flying 
qualities. Additionally, the flight test results provide a ref- 
erence base for the correlation of theoretical prediction and 
wind tunnel experimental results. 

Much of the research effort referenced above was accomplished 
and the results were published in numerous isolated reports, 
during and shortly after the end of World War II. Immediately 
after this period, interest in the problems of straight wing air- 
craft was diverted to the pressing problems of supersonic mili- 
tary aircraft incorporating swept wing configurations. In general, 
unswept (straight) wing technology is not applicable to swept 
wing configurations. Consequently, the wealth of information 
pertaining to the maximum lift and stalling of unswept wing air- 
craft has not previously been coordinated and hence has not re- 
ceived adequate attention. Nevertheless, in the general aviation 
field, interest still centers on .subsonic aircraft incorporating 
unswept wing configurations. A need therefore exists for the 
application of available unswept wing stall technology to the 
design of such aircraft. 

This report presents a comprehensive bibliography of prior 
work in the field, and, insofar as practicable, presents the most 
pertinent information in a form suitable for design application. 

In the preparation of this report a comprehensive review has been 
conducted of all pertinent literature although no attempt is made 
to synopsize each of the reports under one section. Effort has 
been made, however, to incorporate the information gained from 
this extensive review into design guidance procedures, recommen- 
dations, cautions, etc. Based on this comprehensive review of 
pertinent literature a mathematical model was formulated and 
programmed for the CDC 6600 digital computer. The computer which 
employs available nonlinear wing section characteristics can be 
utilized to predict maximum lift and the spanwise load distribu- 
tion of a wing with or without fuselage. 
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As a result of this study it can be concluded that design 
for good airplane stalling characteristics is still part art and 
part science. It appears obvious, however, that application of 
the available knowledge in early design stage will greatly improve 
the probability of obtaining satisfactory stalling characteristics 
or, at least, may yield an airplane design whose characteristics 
can be made acceptable as the result of minor modifications during 
early flight test phase. 



SECTION 2 


BASIC CONSIDERATIONS OF AIRPLANE STALLING 
2.1 MINIMUM REQUIREMENTS FOR ACCEPTABLE STALL CHARACTERISTICS 


The MINIMUM REGULATORY requirements for the stalling behavior 
of small aircraft are stated in the Federal Aviation Regulations, 
Part 23 "Airworthiness Standards: Normal, Utility and Acrobatic 

Category Airplanes" (Reference 10). With the legalistic qualif- 
ications deleted, the regulations require that for specified 
power, gear and flap settings acceptable stalling characteristics 
be demonstrated for two flight maneuvers, one in straight flight 
with the wings level and one in a coordinated turn. In both cases 
the primary control manipulation is a steady progressive upward 
movement of the elevator until the airplane is stalled or the 
elevator reaches its stop. 

The demonstration procedures and acceptable stalling charac- 
teristics are defined as follows: 

For the straight-flight maneuver, the airplane is trimmed 
at a speed fifty percent greater than the stalling speed before 
the elevator movement is started. With the usual three control 
system, it must be possible to produce and correct roll by unrev- 
ersed use of the rolling control and to produce and correct yaw 
by unreversed use of the directional control up to the time when 
stall becomes apparent, i.e. when an uncontrollable downward 
pitching motion develops or until the elevator reaches its stop. 

For the two-control system airplanes, rolling motions must 
be produced and corrected by unreversed lateral control without 
excessive yaw. During recovery from the stall, it must be possi- 
ble to prevent the occurrence of more than 15 degrees of roll or 
yaw by normal use of the controls. 

For the turning-flight maneuver, the airplane is placed in a 
steady, level, coordinated turn with a 30-degree bank angle before 
the speed is reduced. When stall occurs, it must be possible to 
regain normal level flight without excessive loss of altitude 
or uncontrollable rolling or spinning tendencies. 

For both maneuvers as stalls are approached, there must be a 
clear and distinct stall warning which begins at a speed about 5 
to 10 miles per hour higher than the stalling speed and continues 
until the stall is reached. For airplanes that cannot be stalled 
for the two maneuvers specified, it must be shown that if they 
can be stalled in steep climbs, then recovery shall not require 
excessive speeds and/or accelerations. Multi-engine aircraft 
have additional single engine-out requirements pertaining to 
stalls in turning flight. 
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In summary , the regulations require that for gentle maneuvers 
the aircraft shall have some stall -approach warning, effective 
lateral control up to the stall, and sufficiently effective lat- 
eral control after the stall to restrict yawing and rolling dis- 
turbances to small angles. Furthermore, recovery from the stall 
shall not involve excessive altitude loss, speed increase, or 
st rue tur al 1 o ading . 

2 • 2 desirable STA LL CHARACTERISTICS 

For more safety, it would be better if there were no uncon- 
trollable motions during a stall and recovery. A 15-degree change 
in the angle of bank or heading could have serious consequences if 
it occurred close to the ground and off the end of the runway 
following a misjudged landing and the start of the go- around climb. 

Since it is generally conceded that there is no such thing as 
a "good” stall (Reference 11), the safest low-speed characteris- 
tics are those that do not change noticeably to the limit of the 
up-elevator travel. These are, of course, the characteristics 
sought by "stall-proofing" an airplane, such that the angle of 
attack for maximum lift can never be attained. This is not easy 
to accomplish for reasons which will be discussed later. 

For airplanes which can be stalled (i.e. aircraft angle of 
attack exceeds that for maximum lift), there is no unanimity of 
opinion on the least desirable type of motion that may result. The 
consensus, however, is that probably the least undesirable stall 
characteristic would be for the nose of the airplane to drop ab- 
ruptly by a small but distinguishable amount before the occurrence 
of excessive rolling and yawing motions requiring the use of either 
lateralLor directional controls. Also, under such, conditions, if 
the elevator control is eased forward the airplane should promptly 
return to unstalled flight. Slight skidding or yawing motions 
should have no appreciable effect on these characteristics. 

The preceding are the more or less limiting conditions. Any- 
thing between the case of mild nose dropping and the stalling 
characteristics defined by the Federal Aviation Regulations (FAR) 
is a matter of the degree of violence of the rolling and yawing 
motions that develop at the stall. It should be noted that aside 
from the specified angular deviations in yaw and roll, the wording 
of the regulations is qualitative and open to individual inter- 
pretation. 

2 . 3 ST ALL- PROOF ING 


Stall-proofing aircraft by limiting the up-elevator travel 
was suggested at least thirty-five years ago and is a subject 
on which much research time and effort has since been expended. 
There are two major reasons why this solution has not been more 
widely applied. 
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The first, but not necessarily the more important of these 
reasons, is the opinion that a stall-proof airplane may not be 
completely acceptable nor most saleable to the specific segment 
of the general public most interested in personal flying. This 
group is considered to be composed largely of more adventurous 
persons (e.g. the buyers of sports cars, the water skiers) who 
want more than transportation from an airplane. Therefore, the 
amount of effort that might be put into designing an airplane to 
be stall-proof depends on the purpose and the market for which the 
airplane is intended. It is possible to stall-proof an airplane 
for certain flight conditions and not for others. Fortunately 
the landing approach condition is one of the simplest and one for 
which some effort may be warranted. 

The other reason why most aircraft are not completely stall- 
proofed is the technical difficulty of so doing. The elevator 
angle for a given angle of attack varies with, a large number of 
interrelated aircraft dimensional and mass parameters. Among the 
more important not under the control of the designer are the 
center of gravity location on a given flight, the throttle, flap 
and trim settings; and the exterior surface condition, particular- 
ly the wing and horizontal tail surfaces. Added to these are all 
those other parameters among which the designer normally has to 
compromise. The basic longitudinal balance and, hence, the eleva- 
tor angle for stall, depends on the relative proportions of the 
wing and the horizontal tail surfaces, the tail length, vertical 
location of the tail, the orientation of the tail relative to the 
propeller slipstream, and the proximity of the tail to the ground 
during takeoff and landing. 

2.4 STALL WARNING 


It is less difficult to obtain stall warning in the airplane 
as compared to stall proofing. Although it may not always be 
possible to provide for adequate stall warning in the design stage, 
it can in many cases be incorporated after aircraft construction. 

There are a number of airplane characteristics which can 
serve as a means of providing warning to the pilot of his approach 
to the stall. Some of these characteristics are large increases 
in the elevator force and in the control stick travel near the 
stalling speed, initially controllable rolling and pitching motions 
and vibrations or ”buf feting” of the airplane and the control stick 
The variation of control forces which are functions of elevator 
hinge moments are too dependent on the friction and elasticity in 
the control system between the elevator and stick to provide the 
consistency needed for stall warning. On the other hand an appre- 
ciable increase in the elevator control travel as stall is 
approached, airframe buffeting and initially controllable angular 
motions have each been found to give acceptable stall warning. 

The difficulty is to provide for the start of these warnings in 
the required speed range and to assure that the airplane angular 
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motions do not become too severe in amplitude before the stall- 
ing speed is reached. 

Although inherent stall-warnings are to be preferred, stall 
warnings can be synthesized by detecting consistent changes in 
the airflow about the airplane and using the information to 
activate a mechanical or electronic device to alert the pilot of 
his approach to stall. The earliest form of a synthesized stall 
warning was probably the stall red-line on the airspeed indicator. 
The stall airspeed, however, is one of the poorer indicators of 
impending stall because it depends on aircraft weight and conse- 
quently varies with passenger load and fuel consumption. Varia- 
bles which are directly related to the stall such as the airplane 
angle of attack and surface pressure differences give much more 
reliable indications of impending stall than does airspeed. 

Synthetic stall-warning devices can be designed to alert the 
pilot of his approach to the stall through his sense of sight, 
hearing or feel. Neither of the first two have been found very 
satisfactory. In the case of sight, the pilot’s attention during 
takeoff and landing when unintentional stalling is most probable 
is required outside of the cockpit. In the case of hearing, 
there is chance for confusion as horns have long been used as a 
reminder to lower the landing gear. The most satisfactory stall 
warning device would therefore be the one that would alert the 
pilot of his approach to the stall through the sense of feel. 

A description of the most common stall warning devices and their 
principle of operation is presented in Reference 12. 

2.5 FACTORS AFFECTING WING STALL 


The Federal Aviation Regulations define the stall character- 
istics in terms of movement of the elevator and of the effective- 
ness of the aileron and rudder controls. The wing, however, is 
the primary element of the airplane affecting the stalling 
characteristics. Its size relative to weight sets the stalling 
speed, and its proportions determine the location of the stall, 
its rate of propagation, and the violence of the resulting motions 
when the stalling speed is reached. 

In the design of airplanes it is normal to start with the 
wing and select its dimensions by a series of compromises so as 
to assure good performance, acceptable stalling characteristics, 
and low structural weight. While not a technical factor a pleas- 
ing appearance must be considered because of the effect of style 
on the marketability of the design. The airfoil section is 
fundamental to wing design and should be considered early in 
the preliminary stages. The aerodynamic data for the selected 
airfoil sections are then used in analyses to determine the 
span loading of the three dimensional wing. It is on the basis 
of computed span-load distributions that a decision is made on 
the final wing proportions. Sections 4 and 5 of this report 
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present a summary of the theory and a description of the computer 
program respectively for determining the wing span load distrib- 
utions required for preliminary design purposes. 

2.6 AIRFOIL SECTIONS 


The selection of the optimum airfoil sections for a wing is 
generally made after consideration is given to the following: 

a) Section profile drag coefficient should be as low as 
possible over a range of lift coefficients near the required 
cruise lift coefficient. 

b) Section maximum lift coefficient should be as large as 
possible since this has a direct bearing on the maximum value of 
the overall wing lift coefficient which in turn governs the 
stalling speed. 

c) The chosen sections should have sufficient depth to 
contain the wing structural members and other items such as 
landing gear. 

d) Stalling characteristics of sections near the wing tip 
should be gradual so as to avoid the possibility of sharp wing 
drop. 


e) The section critical Mach number should be as high as 
possible to avoid transonic drag rise, if a high speed aircraft 
is being considered. 

All of these requirements cannot be satisfied by any one 
airfoil section and some compromise must be made. For example, 
the need for sufficient structural depth, e.g. large thickness- 
chord ratios, conflicts with the requirement for high values of 
maximum lift coefficient. Again, a section chosen for its 
high value of maximum lift coefficient might have undesirably 
sharp stalling characteristics. 

While section maximum lift coefficient is of particular 
importance in regard to wing stalling characteristics, the 
section having the greatest value of maximum lift might be 
sensitive to small dimensional deviations and hence might not 
realize its ideal performance. 

Wing stall characteristics can also be affected by air 
turbulence or gustiness. The sideslipping and yawing motions 
produced by gustiness can induce flow breakdown for certain 
wing sections thus resulting in large changes in wing lift. 

Much work has been done to correlate flow patterns at the stall 
with the geometric properties of wing sections and the shape 
of the lift curve. There are three representative types of 
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airfoil section stall generally considered in the technical 
literature, e.g. Reference 13. These types are: 


a) 

Type 

1 

Trailing Edge Stall 

b) 

Type 

2 

Leading Edge Stall 

c) 

Type 

3 

- Thin Airfoil Stall 

The 

lift 

curves for the three types (reproduced from Ref 


erence 13) are compared in Figure 1 and a discussion of each 
type of stall is presented below. 

2.6.1 Type 1 - Trailing Edge Stall 

The peak of the lift curve for this type of stall is charac- 
teristically rounded with a maximum value of lift coefficient 
of approximately 1.5. The loss of lift after the stall is 
gradual which is considered to be the least undesirable type 
of stall. The flow for this type of stall is characterized by 
a progressive thickening of the turbulent boundary layer on the 
upper surface of the airfoil as the angle of attack is increased. 
This is followed by an eventual flow separation which starts at 
the trailing edge and gradually moves forward as the airfoil 
angle of attack is increased from about 10o # Maximum lift is 
obtained when the point of separation reaches about the 50 per- 
cent chord station. Beyond maximum lift the forward progression 
of the separation point continues at about the same rate as prior 
to stall. The trailing-edge type of stall is generally associated 
with wings approximately 15% thick or greater. 

2.6.2 Type 2 - Leading Edge Stall 

Airfoils with the leading-edge type of stall usually have a 
slightly greater maximum lift than those with a trailing edge 
stall. At the stall, however, there is a large sudden drop 
in lift associated with an abrupt separation of the flow from 
the upper surface near its leading edge. The separation is 
attributed to the behavior of the upper surface boundary layer. 

At an angle of attack well below that for maximum lift the 
boundary layer, which is laminar at the time, separates from the 
upper surface quite close to the nose, becomes turbulent and 
reattaches to the surface. The localized region of separated 
flow is referred to as a "laminar separation bubble”. Up to 
the angle at which the bubble forms the boundary layer behaves 
the same as for airfoils having trailing-edge separation. The 
difference at higher angles is due primarily to the nose radius. 

As the angle of attack is further increased, the laminar 
separation point moves forward where the local curvature is 
greater. At stall the main flow cannot expand rapidly enough 
for reattachment and a sudden and complete disruption of the 


9 




Figure 1. 


Representative Lift Curves. 

(Reproduced from Reference 13) 



flow takes place. As the nose radii of airfoils normally 
decrease with thickness leading-edge stall is more likely to 
occur with airfoils with thickness ratios of 12% or less. This 
type of stall is considered to be the most undesirable. 

2.6.3 Type 3 - Thin Airfoil Stall 

This type of stall occurs on all sharp leading-edge airfoils 
(regardless of thickness ratio) and may also be encountered on 
rounded leading edge airfoils with sufficiently low thickness 
ratios e.g. t/c<9%. As can be noted from Figure 1, the airfoils 
with Type 3 stall have substantially lower values of maximum lift 
coefficients than those associated with the airfoil types 
discussed above. 

The type 3 stall is characterized by flow separation from 
the leading edge with subsequent reattachment further downstream. 
The exact mechanism of flow reattachment is not clearly understood. 
Experimental observations have shown that at low angles of attack 
the flow reattaches to the upper surface of the airfoil at a 
short distance behind the leading edge and stays attached up to 
the trailing edge without further separation. With increase in 
angle of attack, the point of flow reattachment progressively 
moves toward the trailing edge and the stall is fully developed 
when the reattachment point occurs approximately at the trailing 
edge. 


Figure 2, reproduced from Reference 14, more clearly delin- 
eates the types of stall discussed above and indicates some basis 
for a correlation between stalling characteristics, Reynolds 
number and leading edge shape. 

2.7 WING PLANFORM EFFECTS 


In the design of the wing for stalling characteristics, wing 
planform is as important as the characteristics of the airfoil 
sections. Probably the most important planform parameters are 
the aspect ratio, taper ratio, and the sweep and twist angles if 
either is used. They have an important influence on airplane 
performance and wing weight as well as the stalling characteris- 
tics. Therefore, effective compromises between these parameters 
must be made to obtain the required performance, handling qualities 
and structural integrity of the airplane. 

The aspect ratio for a given wing area will directly affect 
the wing span. Low aspect ratio and short span are conducive 
to maneuverability for an aircraft intended for aerobatics. High 
aspect ratios are attractive from the point of view of 
maximum aerodynamic efficiency or L/D. However, an increase in 
aspect ratio generally results in an increase of wing root 
bending moments and probably wing weight. One fact that should 
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Reynolds number 



Upper-surface ordinate at O.OI25 chord, percent chord 


Figure 2. - The Low- Speed Stalling Characteristics of 
Airfoil Sections Correlated With Reynolds 
Number and the Upper- Surface Ordinates of 
the Airfoil Sections at the 0. 0125-Ch.ord 
Station. (Reproduced from Reference 14). 
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be considered relative to the compromise between maximum L/D 
and wing root bending moment is that the lift coefficient at 
which maximum L/D is attained also increases with aspect ratio. 

The performance characteristic most affected by maximum L/D is 
the aircraft range. For personal aircraft with relatively 
light wing loadings operating at altitudes below 10,000 feet, 
the speed for maximum range may be too slow to make it an 
attractive cruise speed. The increase in the value of the 
maximum L/D obtained by increasing the aspect ratio, therefore, 
may be of little practical use. In fact, if the L/D is much 
increased for the lift coefficients used in the approach glide, 
the effect may be to flatten the approach path, make judgment 
of the landing point more difficult and give the airplane a 
tendency to "float" after the landing flare. This statement is 
intended only to imply that increases in aspect ratio do not 
necessarily lead to improved aerodynamic performance. Improved 
performance may be better obtained by a reduction of fuselage, 
engine and landing gear drag. 

Taper is normally employed to increase the wing chord and 
spar depths at the wing root and, hence, offset the adverse 
weight effects of increasing the aspect ratio. At the same time 
it increases wing torsional rigidity which is an important 
structural consideration. Wings can be tapered either in 
thickness or planform or in any combination of both. When com- 
bined taper is incorporated in the wing and straight line surface 
fairing is used for structural simplicity, airfoil sections be- 
tween the fairing stations may be slightly distorted. 

Sweep shifts the aerodynamic center of the wing fore-and- 
aft relative to the wing root. Small amounts of it can therefore 
be used to adjust the longitudinal aerodynamic balance or to 
improve the position of the aircraft center of gravity relative 
to the main landing gear. 

Wing twist affects the span load distribution by changing the 
local section angle of attack along the wing span. 

Negative wing twist or washout may be used effectively to move 
the stalling point inboard. 

High-lift devices such as leading edge slots or trailing 
edge flaps, if they are installed across the full span, influence 
the wing stall only through their effect on the section lift 
characteristics. Rarely, however, are either leading edge slots 
or trailing edge flaps used over the entire wing span. Leading 
edge slots over the tip portion of the wing were one of the first 
means tried for improving the wing stalling characteristics. Such 
slots increase the angle of stall of the tip sections of the 
order of ten degrees. Much experimentation with relatively 
simple fixed slots, as well as with the movable Handley-Page type, 
led to the early conclusion that neither the complexity 
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and cost of the automatic slots nor the drag penalty of the 
fixed slots were warranted for small personal aircraft. There 
has been a recent revival of interest in movable leading edge 
devices for the tip portions of the wings of high-speed, swept- 
wing commercial and military aircraft where their use can be 
economically justified. 

Flaps are usually installed inboard of the ailerons. They 
depress the angle of attack for zero lift of the airfoils to 
which they are attached while having only a secondary effect on 
the stall angle. Inboard flaps, therefore, have the effect of 
a discontinuous wing twist with wash-in starting at the flap end. 

Application of the simple lifting line theory as described 
in Reference 15 to predict the stalling characteristics of 
flapped wings will indicate a discontinuity in the spanwise 
distribution of maximum section lift coefficient and an assoc- 
iated initial stall at the end of the flap. Inasmuch as this 
discontinuity does not exist in three dimensional flow, early 
attempts to predict the stalling characteristics of flapped 
wings (Reference 16 and 17) by the method of Reference 15, did 
not yield good correlation with experiment. Subsequently, 
a far more reliable method of applying the lifting line theory 
was developed (Reference 2). The latter method predicts the 
point of initial stall quite well and is described in detail 
in Section 4 of this report. 

2.8 CONSTRUCTION TOLERANCES AND SURFACE IRREGULARITIES 


In most cases the final wing design represents a compromise 
among the best dimensions for performance, stall character- 
istics and structural weight, although occasionally some 
importance may be given to appearance. Curved leading edges, 
taper and sweep can add to the attractiveness of a design and 
give the impression of speed. While such features so used may 
not necessarily penalize the design, their impact on flight 
safety and construction costs should be carefully studied. 

No adequate method of defining either the permissable 
deviations of the airfoil ordinates from their specified values 
or the permissable degree of surface roughness and waviness has 
been devised. Reference 5 gives some general guide lines 
concerning the surface conditions and construction tolerances 
that must be adhered to in order to achieve extensive regions 
of laminar flow at low values of the lift coefficient. In 
general, it would appear that if values of c / max comparable to 
those presented in Reference 5 are to be achieved, the surface 
conditions and construction tolerances over the forward upper 
surface of the airfoil (approx. 0 to 10% chord) should be at 
least as good as those indicated in the reference. As shown 
by the results of Reference 18 "very small errors in airfoil 
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contour, particularly around the leading edge, could cause 
large changes in the stalling angle of attack and the resulting 
value of maximum lift coefficient." 

The advisability of applying the stall analysis procedure 
in the early stages of the airplane design cannot be over- 
emphasized. This stall analysis procedure relies on the use 
of published airfoil section data characteristics. The designer 
is warned, therefore, that unless the aircraft as finally 
constructed incorporates airfoil sections that are approximately 
equal in contour and surface condition to those utilized in 
obtaining the two-dimensional airfoil section data, the calculated 
stalling characteristics will not necessarily be comparable to 
those encountered on the flight vehicle. 

2.9 WING -FUSELAGE FAIRINGS 


One of the more important structural considerations affecting 
aircraft stall is the fairing between the wing and fuselage. 

The method of Multhopp used for calculating the effect of the 
fuselage on the wing lift only accounts for changes in local flow 
angles of the root sections due to presence of the fuselage. 
Neither this method nor any other method considered in the course 
of this study adequately treats the modification of the surface 
pressures on the wing and fuselage as a result of the mutual 
interference. Such changes seriously influence the boundary layer 
at the junction and can cause premature separation of the main 
flow. The turbulence in the wake behind the separated region 
may intersect the tail plane where it may induce vibratory loads 
which if mild, can be used as stall warning. It may, however, 
result in vibrations severe enough to impair the structural 
integrity of the aircraft. The violence of the vibrations depends 
on the extent of the separated region and the location of the 
horizontal tail relative to wake. 

Wing fuselage flow interaction effects largely depend on the 
position of the wing relative to the fuselage. For the high 
wing position, where the influence of the fuselage is confined 
to the lower less critical surface of the wing these effects 
are quite small. The low wing position with the wing tangent to 
the fuselage introduces the largest wing-fuselage interaction 
effects which increase as the angle between the wing-fuselage 
surfaces becomes more acute. 

The wing-fuselage flow interaction can be to some degree 
controlled by the shape and size of the root fairings associated 
with each wing position. The high wing position dictates a 
minimum amount of fairing from aerodynamic considerations, 
whereas for the low wing position a considerable amount of 
fairing may be required to ensure proper flow conditions on the 
upper surface of the wing close to the fuselage. Also, careful 
consideration should be given to the chord-wise shaping of 
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the wing root fairing, since it effects the expansion of the 
flow over the after portion of the junction. The chord-wise 
shape of the root fairing largely depends on the fore-and-aft 
position of the wing relative to the fuselage maximum diameter. 

Since the geometry of the wing-root fairings critically 
depends on detail design of the junction there are no established 
theoretical methods by which the size and shape of the wing- 
root fairings can be predicted. Some empirical design rules, 
which can be employed for this purpose are presented in 
References 19. However, it is generally recommended that 
the final size and shape of the wing-root fairings be determined 
by an experiment either in a wind tunnel or in flight. 

The fairing problems discussed above are similar to those 
associated with the nacelles of wing -mounted engines. Location 
of the nacelles close to the fuselage may result in some 
additional interaction effects of the nacelle-fuselage flow 
fields. However, the propeller radius requirements generally 
ensure sufficient spacing between the nacelles and the fuselage 
to prevent serious interaction. 

2.10 PROPELLER SLIPSTREAM CONSIDERATIONS 


Another important consideration of airplane stall is the 
effect of the propeller slipstream in power-on flight. This 
effect is introduced through an increase of the local velocity 
over the wing immersed in the slipstream and the change of wing 
local angle of attack due to slipstream rotation. 

The increased velocity tends to stabilize the flow over the 
wing immersed in the slipstream. The rotation within the slip- 
stream tends to increase the local angle of attack of the wing 
sections behind the upgoing propeller blades and decrease the 
local angle of attack behind the downgoing blades. The overall 
effect is usually that of promoting an asymmetrical stall. The 
stall of that portion of the wing behind the upgoing blades is 
hastened whereas that behind the dowi going blades is delayed. 

In the case of a single engine airplane the asymmetrical stall 
can promote serious wing dropping tendencies when the airplane 
is operating in the vicinity of CLmax- Reference 8 presents the 
results of extensive experimental observations of the effect of 
propeller operation on wing stalling. It is interesting to 
note, from the results presented therein, that in the case of 
one of the single engine airplanes investigated, the action of 
the propeller in promoting an asymmetrical stall is more adverse 
at the engine idle power condition (Tc = 0) than at the power-on 
condition investigated (Tc = 0.2). In the case of another single 
engine airplane investigated, the effect of power was in the 
reverse order. 



The effects of the asymmetrical stalling in the case of 
multi-engine aircraft can not be clearly defined primarily due 
to the fact that these effects are largely dependent on a spec- 
ific combination of various geometric parameters related to each 
aircraft. 

2.11 STABILITY AND CONTROL CONSIDERATIONS 


In the above discussion of airplane stall characteristics, 
it is assumed that the aircraft is in a steady unyawed trimmed 
flight, and that the stall is developed through a gradual increase 
of airplane angle of attack. However, the steadiness of actual 
flight depends on the longitudinal stability of the aircraft, the 
control effectiveness, the rate at which the pilot moves the 
control stick and the atmospheric turbulence. 

With reference to atmospheric turbulence, little can be done 
in design to reduce its effect on the stall. Even assuming 
that the airplane stability and control are adequate, severe tur- 
bulence can upset any aircraft and therefore it should be avoided 
where possible. In an unavoidable flight in turbulent weather, 
the Federal Regulations require the pilots of commercial airlines 
to fly the specific airplane at its designated speed low enough 
to minimize the structural loads due to the gusts, but high 
enough to assure that the angle of attack change produced by the 
probable maximum upward gust will not stall the aircraft. The 
designer of the private airplane can only suggest that his cust- 
omers follow a similar practice. 

In addition to the angle of attack changes, gustiness will 
induce rolling and yawing motions that can result in premature 
unsymmetrical wing stall and sudden roll-off. The danger depends 
on the stall margin of the wing and on the lateral stability of 
the aircraft. One of the requirements for good lateral stability 
is an adequate fin area. The relation between the fin area and 
the wing dihedral angle determines whether the yawing or rolling 
component of the lateral oscillation predominates. The occurrence 
of light oscillations as stall is approached has been suggested 
as an acceptable stall warning. This phenomenon is usually 
associated with a reduction of fin effectiveness at large wing 
angles of attack and is caused by blanketing of the lower portion 
of the fin by the fuselage or submersion of the fin in a low- 
energy wing wake . 

The longitudinal stability and the elevator effectiveness 
combine to determine the stick movement required at different 
parts of the speed range. In general, high-wing positions result 
in more longitudinal stability at low speeds than at high speeds. 
Low-wing positions produce the opposite effect. The elevator 
effectiveness varies with the relative proportions of the eleva- 
tor and stabilizer, and it decreases after the angle of the 
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elevator exceeds a critical angle. Hence, for an increasing rate 
of stick movement as the stall speed is approached, as suggested 
for stall warning, the high-wing configuration with a small 
elevator is more desirable. 

There are, however, a considerable number of other factors 
affecting the stability and control characteristics which in turn 
influence the airplane stall. The subject is too complex for a 
short comprehensive treatment in this report. Among others, the 
aircraft stability and control is affected by the size, location 
and type of flaps; the wing wake and its position relative to 
the horizontal tail; the orientation of the propeller slipstream 
relative to both the horizontal and vertical tails and the vert- 
ical, as well as the fore-and-aft, position of the center of 
gravity. These and other factors are adequately treated in the 
available technical literature, e.g. References 20 and 21. 

2.12 FLIGHT VERIFICATION AND CERTIFICATION 


The theoretical methods and the design procedures, presented 
in this report should reasonably well predict the stalling 
characteristics of unswept wing aircraft, particularly for gliding 
flight such as in the landing approach. Since these methods and 
procedures are based upon a number of simplifying assumptions, 
it is always advisable to verify the theoretical predictions by 
actual flight tests. 

The test procedures to be followed are specified in the 
Federal Aviation Regulations for Certification Demonstrations 
cited earlier. If the tests show a need to improve the aircraft 
stalling characteristics, a flow visualization technique, utiliz- 
ing tufts is often employed in gaining an understanding of the 
flow deficiencies. The tufts which are normally short lengths 
of string attached at intervals on the upper surface of the wing 
and sometimes on the sides of the fuselage, can be observed and/ 
or photographed from within the test aircraft or the chase plane. 
The flow patterns indicated by the tufts will define areas of 
attached and separated flows and can show the spanwise or chord- 
wise location of initial flow breakdown and the rate at which 
separation increases with small decreases in speed. 

The tuft studies, the motions after initial stall and the 
stall-related vibrations felt in the structure and the control 
system will generally indicate the design changes required for 
improving the stalling characteristics after the airplane has 
been built. No attempt will be made to list or discuss all the 
possibilities in this report. 

One possibility is that a part of the wing will stall sudden- 
ly with almost no warning. This type of stall can occur with 
airfoils exhibiting leading edge stall and with wing proportions 
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yielding low stall margins at the sections close to the point 
of initial flow breakdown. Interference between the flows over 
adjacent sections promotes separation which can spread across 
the wing very rapidly. The means most used for improving the 
stalling characteristics in such a case is the installation of 
a triangular- shaped projection or spoiler, Figure 3 (a), for a 
short distance along that part of the leading edge over which an 
earlier stall is desired. Attainment of satisfactory stalling 
will generally require a trial- and-error procedure, varying the 
span and spanwise location of the spoiler and its radial loca- 
tion at the leading edge. It may be difficult to obtain a con- 
figuration of the spoiler which is effective with full power but 
not too effective when power is reduced. Spoilers of this type 
can increase the power-off minimum speed by an objectionable 
amount and produce large vibrations of the wing, stabilizer or 
fuselage. 

The other extreme possibility is a root stall starting well 
above the minimum speed. The associated thick turbulent wake, 
if it impinges on the tail, will produce stall warning, but may 
also make it impractical to fly the airplane in the speed range 
required for the landing approach. If the root fillet is prop- 
erly designed it may be necessary to modify the wing leading 
edge of the root sections to increase their forward camber and 
reduce their effective angles of attack. One such, leading edge 
modification is shown in Figure 3(b). It is frequently necessary 
to tailor the leading edge modification such, that a small vibra- 
tion remains for stall warning. Such modification may be util- 
ized to control the asymmetrical stall that sometimes occurs 
within the propeller slipstream. 

The actual cases are rarely as clearcut as the two just 
discussed and the reason for poor stalling characteristics is 
often not clearly understood. During World War II, many aircraft 
with configurations similar to those of current private owner 
designs had to undergo post-design stall improvement testing. 

Many of these cases are in the literature and have been included 
in the bibliography presented in this report. 
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(a) Sharp Leading Edge Strip to Hasten Stall. 



Cb) Modified Nose Radius and Camber to Delay Stall 


Figure 3. - Wing Leading Edge Modifications 
for Controlling Wing Stall. 
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SECTION 3 


THEORETICAL ANALYSIS 


Presented in this section is a brief review of the available 
theoretical analysis which formed the basis for developing the 
mathematical model and the computer program contained in this 
report. Specifically, in selecting the most suitable theoretical 
approaches, due consideration was given to the past work in the 
fields of wing theory and wing-body interference theory. The 
significant contribution in these areas and their applicability to 
the present program are discussed in the following pages. 

3.1 REVIEW OF THE AVAILABLE THEORIES 


3.1.1 Wing Theory 

The simplest wing theory involves the concept of the lifting 
line whereby the wing is replaced by a vortex filament fixed in 
the wing and a system of trailing vortices. The vortex filament 
is known as the bound vortex or lifting line. At each point on 
the span the strength of the vortex is proportional to the local 
intensity of lift and is related to it through the Kutta-Joukowski 
law, which is based on the requirement of smooth flow from the 
wing trailing-edge. In accordance with the Helmholtz theorem of 
vortex continuity it is assumed that with each elemental change 
in the spanwise distribution of the strength of the bound vortices 
there is associated a free vortex which is shed at the wing trail- 
ing edge and which passes downstream with the general mass of the 
fluid. The strength of this vortex is equal to the incremental 
change in spanwise circulation. This system of free or "trailing” 
vortices induces velocity components (called downwash) at the wing 
and thus causes a change in local angles of attack of each wing 
section. 

The problem of evaluating the downwash at each point is 
difficult because of the interrelation of downwash, lift distrib- 
ution and wing planform. Using the lifting line theory Prandtl 
(Reference 22) obtained a solution for the case of a wing with an 
elliptical lift distribution. Glauert (Reference 23), using 
Fourier analysis, developed methods for obtaining solutions for 
wings of any planform and twist. These methods were used by 
Anderson (Reference 15) to determine the characteristics of wings 
for a wide range of aspect ratios, taper ratios and linear twist 
distributions . 

All of these approaches discussed above involve the assumption 
of a linear variation of section lift with angle of attack. This 
assumption was not utilized in the iterative methods developed 
by Sherman, Tani, Multhopp, and Boshar (References 24 through 27) 
which employ nonlinear section lift data in the computations. 
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Based on these iterative methods Sivells and Neely 
(Reference 1) developed solutions which yield excellent agree- 
ment with the test data up to angles-of -attack close to stall. In 
a later paper Sivells and Wes trick (Reference 2) extended the 
method of successive approximations to the calculation of the 
aerodynamic characteristics of wings with deflected flaps or 
ailerons. 

All of the above solutions and methods which are essentially 
based upon the Prandtl lifting line theory are known to be 
inadequate for wings of aspect ratio less than about 3. For wings 
of this class the influence of the chordwise loading can no 
longer be neglected and resort must be made to the general theory 
of lifting surfaces to obtain solutions. 

Lifting surface theory involves finding a potential flow 
solution which satisfies the Kutta condition all along the span 
while at the same time satisfying the boundary condition that 
there is no flow through the wing surface. 

Solutions of varying complexity and accuracy have been 
advanced by many authors. Early attempts were made by Weissinger, 
Mutterperl, and Schlichting (References 28, 29 and 30 to use 
simplifying physical models for the approach to the general 
problem, e.g. placing a lifting line at the quarter-chord point 
and satisfying the downwash condition at the three-quarter-chord 
position. Falkner (Reference 31) proposed a vortex-lattice 
treatment of the wing thereby approaching a truly continuous 
lifting surface. 

Attempts to use a continuous lifting surface theory without 
resort to arbitrary physical assumptions or models are 
exemplified by the work of Garner (Reference 32 and 33) and most 
notably by that of Multhopp (Reference 34). 

When lifting surface theory is used to predict load distrib- 
utions on high aspect ratio unswept wings at low angles of attack, 
the results do not differ significantly from those computed 
using lifting-line theory. 

However, to date, lifting surface theory has not been 
successfully modified to permit the use of nonlinear section lift 
data and hence cannot be expected to give reliable predictions 
of load distribution at wing angles-of -attack near the stall. 

For this reason the lifting-line theory, as modified and 
presented in Reference 2 has been chosen as the more 
appropriate method and the one which is better suited to the 
present task. 
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3.1.2 Wing-Body Interference Theory 

Many methods exist for the calculation of air loads on 
wing-fuselage combinations and these have been summarized by 
Schlicting (Reference 35), and by Flax and Lawrence (Reference 

36) . In view of the foregoing selection of lifting-line theory 
to calculate essentially wing characteristics, the following 
discussion will be limited to those wing-body interference methods 
which have been developed for use with the lifting-line approach. 

The spanwise lift distribution over a wing-body combination 
of minimum induced drag was first treated by Lennertz (Reference 

37) for a body consisting of a infinitely long circular cylinder. 
His solution was generalized by Pepper (Reference 38) to include 
bodies of any cross-section . 

Multhopp (Reference 4), using a conformal mapping technique, 
obtained a solution for the case of a high aspect ratio wing 
mounted on an infinitely long cylinder of any shape. This method 
may be applied if a function can be found which maps the body 
cross-section conformally onto a circle or a straight line. The 
advantage of the method is that the transformed cross-section 
shape automatically becomes a streamline, thus satisfying the 
fuselage boundary conditions. 

An alternative approach is that first used by Lennertz (Ref- 
erence 37) in which singularities in the form of image vortices 
are introduced within the fuselage cross-section to satisfy, 
approximately, the boundary conditions. Zlotnick and Robinson 
(Reference 39) applied this method to the case of swept wing -body 
combinations with centrally placed wings. 

The use of the image vortex method in combination with the 
lifting-line method of Reference 1 is restricted to wings which 
are centrally mounted on circular fuselages. For wings not 
centrally mounted the equations become extremely complex due to 
the fact that the image vortices no longer lie on the extended 
plane of the wing within the fuselage. For other than circular 
cross-sections, the determination of the number and location of 
the image vortices is difficult. 

Weber, Kirby and Kettle (Reference 40) have modified 
Multhopp r s approach to account for non-zero wing thickness and 
applied this method to low-aspect-ratio swept wings. In view 
of this success and because of its simplicity the Multhopp 
formulation of the wing-body solution (Reference 4) has been 
adopted herein rather than the more complex method of images 
of Reference 37 . 
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3.1,3 Propeller Slipstream Effects 


As pointed out in Section 2.10, the propeller slipstream 
exerts an important influence on wing load distribution which 
in turn affects the aircraft stall characteristics. 

A review of the technical literature indicates that there 
are no adequate theoretical or semi-empirical methods which 
can properly account for the effects of propeller slipstream on 
the spanwise load distribution of the entire wing. Most of the 
available theoretical and experimental investigations, e.g. 
References 41 through 47 predict wing load distributions solely 
for that portion of the wing immersed in the propeller slip- 
stream. This approach is considered to be inadequate for 
predicting stall characteristics of wings spanning the propeller 
slipstream because the portions of the wing outside the 
slipstream cylinder are strongly influenced by the slip- 
stream flow. On the other hand, development of a mathematical 
model for predicting the effects of propeller slipstream on 
wing load distribution within and outside the slipstream 
cylinder is rather complex and is considered to be outside of 
the scope of the present work. 

In view of the above difficulties, it was decided to 
temporarily exclude the power effect from wing stall analysis. 
However, provisions have been made to include such effects at a 
later date, when reliable theoretical methods for predicting the 
effects of propeller slipstream on the load distribution of the 
entire wing become available. Therefore, in developing the 
theory presented herein, the method of Multhopp (Reference 4) 
for wing-body interference effect and the method of Siveils 
(Reference 2) for the isolated wing have been successfully 
combined and are presented in the following analysis. 
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3.2 FORMULATION OF THE ANALYSIS 


3.2.1 Conformal Transformation of Wing-Fuselage Combination 

A wing in the presence of a lifting fuselage (a B >o ) is 
subjected to an upwash which decreases towards the wing tips. 

This upwash has the effect of a variable wing twist, which is 
a function of fuselage and wing geometries, fuselage angle of 
attack and aircraft forward speed. 

The problem of fuselage-wing flow interaction has been the 
subject of numerous investigations in the past. As mentioned 
previously, the approach selected in the present analysis is that 
based on Multhopp’s formulation. In this approach, the wing- 
fuselage combination is conformally transformed into an equiv- 
alent wing with a vertical slit representing the fuselage. The 
vertical slit is aligned with the cross-flow and is therefore 
automatically a streamline. The transformation is applicable 
to fuselage cross-sections which are approximately elliptical 
or which can be transformed to such by the use of other methods. 

If the wing is not centrally mounted on the fuselage it 
transforms into a slightly curved trace whose curvature increases 
with an increase in the vertical distance of the wing relative 
to the center of the fuselage. This introduces an effective 
dihedral to the transformed wing, but its effect is considered 
to be small and is therefore not treated in this analysis. 

Figure 4 shows the geometry of a wing-fuselage combination 
in both the real (u) and transformed (u) planes. The definition 
of real and transformed parameters is presented in the list of 
symbols. 

The wing load distribution, c % (y) , in the u-plane is a 
function of the local effective angle of attack ci e (y) of the 
section under consideration. This effective section angle of 
attack can be expressed as follows: 


c* e (y) = q b + a R + € (y) + A a (y) - a/ (y) 


( 1 ) 


Where a B ,a r , and e (y) are the body angle of attack, wing-root 
incidence relative to the body and the wing section geometric 
twist , respectively. 

The angle due to body upwash A a(y) and that induced by the 
wing trailing vortex system a/ (y) can be determined from the 
following analysis. 
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Figure 4. Definition of Parameters for Transformation 
of Wing-Body Combination 
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Using Reference 4, the wing angle of attack due to body 
upwash for zero wing thickness is given by: 


A a (y) = a B ( R- 57 - 1 ) 


( 2 ) 


where R du/du , is the real part of the derivative of the conformal 
function U (u) . For an elliptical fuselage, this derivative 
can be expressed as follows? 


o 1 

r - 1 f A i _ B 1 V qi2 - e 12 1 

- d u A'-B 1 L , + ei 2 y 2 J 

(a'2- e '2) 


where , 


a' = t (Vy2 + (h-e ')2 +V y 2 +( h +e') 2 } 


(3) 


(4) 


For a circular fuselage equation (3) becomes: 

ho , , A l2 ( y 2 - h 2 ) (5) 

- du (y2 + h2)2 

For wings of non-zero thickness, Reference 40 suggests that the 
upwash at the wing is decreased as compared to that for wings 
of zero thickness. This can be approximately accounted for by 
reducing the wing angle of attack due to body upwash, A a (y) , 
by a factor T taken as constant across the wing span. The 
factor can be expressed as the ratio of the body cross-sectional 
area above and below the wing to the total frontal area of the 
body. Thus, if the wing is not too thick, there follows: 

2 tr yo (6) 

m - 1 V ' 

1 TT A 1 B 1 

Thus, for wings of non-zero thickness, equation (2) can be 
rewritten as: 


A a (y) = T a B ( R-^J- 1 ) 


(7) 


If, however, the conformal function for the thick wing,(R^k-) T 
could be determined, then equation (7) would be Qu 
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A Q (y) - q b 


( 8 ) 


Comparing equation (7) and (8) there results: 




(9) 


Equation (9) represents an approximation to the real part of 
the derivative of the conformal function applicable to wings 
of non-zero thickness. This equation is therefore used to 
relate the induced angle of attack at a point y in the u -plane 
to that at the corresponding point y in the u -plane, thus: 

a i (y) = ®/(y>[e$f] T = a/W [i + T(R-^--d] (10) 


Finally, substituting equations (7) and (10) into equation (1) 
the effective angle of attack for wings of non-zero thickness 
is given by 


a e (y) = a R + <= 



-a;(y)][l + T(R^--D_ 


( 11 ) 


The_only unknown which remains to be determined in equation (11) 
is a/ (y ) . 


The induced angle of attack, in degrees, at a point y , 
the transformed span is given by the familiar relation 


on 


ai = 180 b 


8 7T2 


p b/2 

•H5/2 


,, 5(y)c 
d( B \ 

dy_ 

yi - y 


dy 


( 12 ) 


where c \ (7) is the section lift coefficient in the transformed 
plane at a point 7. It should be noted that since the trans- 
formation is conformal, the circulation T (7) about any wing 
section and the associated section lift coefficient c / (y) in 
the transformed plane are equal to the corresponding values in 
the real plane. Since the geometric quantities of chord length 
(c) and wing twist at each point in the real plane are also 
the same at the corresponding points in the transformed plane 
then 
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Cj( y) = c z (y) = C; 


(13) 


For wings with undeflected flaps, the spanwise lift 
distribution in the transformed plane can be expressed as an 
infinite trigonometric series as follows: 

J ^ = X A n sin(n8) (14) 

b h=l 

where , 

0 = cos-'(^) (15) 


Using equations (12) and (14) there follows, according to 
Reference 23 


a l 


180 
4 7T sin 


00 

£nA n sin(n Q ) 
n=l 


(16) 


Thus, in order to determine the induced angles of attack the 
values of the coefficients A n are required. These are obtained 
as follows: 


The spanwise lift distribution is approximated by a finite 
trigonometric series of r-1 terms corresponding to determining 
the load at an odd number, r-l , of points on the span. That 
is, the range o<Q<ir is divided up into r equally spaced intervals 
where 6 =mir/r m= I , 2, . . . . , r-i . 

Thus, equation (14) can be rewritten as: 


( c i_ c 


m 


r-i 

= X A n sin n 
nr [ 


m 7T 
r 


(17) 


Using harmonic analysis, equation (17) yields the coefficients 
A n of the trigonometric series, thus 



m=i 


c/c 
b m 


sin n 


m 7T 
r 


(18) 


If equation (18) is now combined with equation (16), an expression 
is obtained for the induced angle of attack, thus 
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a» 


r-i 

= _L 80 _/ y 

47Tsing \ 



sin n 


m 7 r 

r 


180 

47 rrsin 0 r 


Z (-T-) 

m=i m n=i L 



Since the induced angle of attack is to be determined at the 
points S at which the load distribution is required, i.e. at 
points q = then 


a ' k 



cos n ( .HP )7 L_ cos n 


or, 


Q/ k = 2 ^ mk 

m-i 


where 


180 


r-i 


m k 4 7T r sin 


Z " 


(k-m)7r (k+m)7r 

cos n cos n 


It can be shown that if k+m is odd then 


/^mk 


180 


47rr sin 


k7T 


I 


I -cos iktmffi 


l-cos 


(k-m) 7 T 


and if k=m 


/^mk 


180 r 


8 7 T sin T- 


finally, when k+m is even and k+m 


( 19 ) 


the 


( 20 ) 


( 21 ) 


( 22 ) 


(23a) 


(23b) 
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/3mk = 0 


(23c) 


The values of /^mk depend only on the number of spanwise stations 
used and are independent of wing aspect ratio and taper ratio. 
Knowing the values of /3 m k * equation (21) can now be utilized 
to determine the induced angle of attack, a/k , at a point on 
the transformed span, in terms of a known or assumed spanwise 
lift distribution. 

3.2.2 Spanwise Load Distribution for a Wing With No Flap 
or a Full Span Deflected Flap 

The method of determining the lift distribution is one of 
successive approximations. For a given body angle of attack a 
distribution of c* c/ b is assumed and the induced angles of attack 
are computed using equation (21). Using equation (11) the 
effective section angles of attack are calculated and the corre- 
sponding values of lift coefficient obtained from airfoil data 
at the appropriate values of the section parameters. This 
process is repeated until the guessed values agree with the 
computed values to within the required tolerance. 

In order to minimize the number of iterations required to 
converge on the final distribution of lift coefficient a system- 
atic method is required to generate increasingly better approx- 
imations. Such a method is developed in Reference 2 and is 
herein presented in slightly different notation. 

For a wing-body combination with zero thickness wing and 
no flap deflection, the basic equation to be satisfied during 
any cycle of the iterative process is the simplified version 
of equation (11), i.e. 


a e (y) = a g (y) + e(y) + a B (R-g)-a/ (y) R-jJ- 


(24a) 


= aq (y)+ e (y)+ a B (R-^- 



or 


a e (y) = a k 



r-i 


t m 


m = i 


b 'm 


^mk 


(24b) 
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a k = a gM + € (y) + a B ( R ) 


(25) 


where 


Assuming linear section lift-curves the relationship for one 
cycle of the iteration is 



-s£ | (-Y-1 » *(ttX * 


(26) 


where is the amount to be added to the approximate value 
ci c/b to obtain the calculated value. 

For the following cycle a value A'k is chosen such that the 
calculated values are equal to the guessed values thus: 




Subtracting equation (27) from equation (26) yields 


(nr) k ( 8 -f-)k I P mk = A k - A' k ( 28 ) 

m=i 

then with 

“'—(nrX*- ( 8 Sk 

= 1 + (~r)k ^ mk (- dir)k m=k 


there follows 


Z G mk A 'm - Ak 
m=i 

j^rnkj 


(29a) 


(29b) 
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where G m |( is the transpose of G^ . 

Thus, the values to be added to one set of approximate 
values to obtain a better approximation are given by: 


{A'm} = [5mk]''m (30.) 

or, 

{ A ' j H K,i ] { Aj } (30b) 


The matrix of coefficients Kn is easily obtained by the 
usual methods of matrix algebra. 

3.2.3 Spanwise Load Distribution for a Wing With a 
Deflected Part -Span Flap 

The deflection of a part- span flap causes a discontinuity 
in the distribution of absolute angle of attack at the end of 
the flap. In order to maintain a continuous spanwise distrib- 
ution of lift a corresponding discontinuity must exist in the 
distribution of induced angle of attack. 

From Reference 2, the complete lift distribution may be 
expressed as the sum of two distributions, thus: 

c '_ c = C M C + § / °12 C \ (31) 

b b ' b ' 

where S(cz 2 C/b) is the distribution due to a unit discontinuity 
in induced angle of attack and c^c/5 is the remainder of the 
lift distribution. Since the latter distribution is continuous 
the multipliers /?mk may be used directly to obtain the correspond- 
ing induced angles of attack. 


a 


* i k 



(32) 


Since ot/ 2k = 8 over the flapped span and a /2k = 0 over the 
unflapped span, the total induced angle of attack is 
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(33) 


at k = a< i k + a< 2k 


Now, if the multipliers were used with the total lift distrib- 
ution 



-)m ^ mk= Z(-^-) m^ mk+ Z S (- 


m=i 


m= i 




mk 


and az 2k were added to both sides, the result would be: 


I -Sm k**2 k = Z (4f) m /3mk > * 2| /f 8 ^-)/S m k 

m=i m = i ° m=i 


(34) 



Rearranging equation (34) yields: 


a h 




(35) 


A comparison of equation (35) with equation (21) shows that 
for a wing with a deflected flap an additional term is required 
which is proportional to the magnitude of the discontinuity. 

This terms represents a correction factor to account for the 
inability of a limited trigonometric series to represent ade- 
quately the spanwise lift distribution of a wing with deflected 
flaps. 


Equation (35) may be rewritten as: 


a ‘ k - 



(36) 
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where the uncorrected induced angle of attack is expressed as 


r-i 


■..■ZC-VV"* 


( 37 ) 


and the correction factor per unit discontinuity is 


a|k = ai|k._y (3 mk 

8 8 4| b 8 


(38) 


The distribution of ( c /2 c /bS) may be expressed in the form: 


c/ 2 c 

bg = 


i 

90 




(cos 9 - cos 9*) lo 9 


r i-cos (0 + 0*) ~ l 7r5*in 9 
1 I -cos (Q-Q*) J 90 




for =1, o < 9 < 9* ; 

for 2 =0, 9*<9< l80 °i 


(39) 


where 9 * = cos 1 (2y^ b) and 9 * and 8 are in degrees. 

The distribution of c /2 c /bS as given by equation (39) depends 
solely on the spanwise position of the discontinuity at the_ end 
of the flap, and applies to outboard flaps which begin at 2y/b = 1 
and end at 27 */b . For a_wing with symmetrical inboard flaps 

extending between 2y/b = 2y*/b and 2y/b =-2_y*/b the distribution 
is obtained by _subtr acting values of c ^2 c /b 8 for 2 y*/b from 
those for - 2y*/b # The correction factor per unit discontinuity 

Gck/8 , given by equation (38) is therefore a function of only 
the spanwise position of the discontinuity and /3 m k . 


It should be noted that two values__ of_ a c/ 8 exist at the 
end_of_the flap, one corresponding to 2y*/b + 0 y (the flap side 

of 2 y */ b ) and the other corresponding to 2y*/b - 0 , (the unflapped 
side). The values are related by: 

(~ g ~t~) 2 y* = * _1 (40) 

b ° b + ° 
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Either of the two values of ac/8 may be used so long as the 
value is used with the proper section lift curve. 

3.2.4 Modification of the Two-Dimensional Section Lift Data 


If the two-dimensional section data were used directly 
discontinuities would be found in the maximum lift distribution 
at the ends of the deflected flaps. Since the discontinuities 
do not actually occur, the two-dimensional data must be adjusted 
to reflect the fact that sections in the vicinity of the end 
of the flap exhibit values of c/ max which lie between those for 
flapped and unflapped sections. Therefore, following the method 
of Reference 2 the two-dimensional data is modified to obtain 
so called three-dimensional section data. 

The root section of a wing with a partial- span flaps deflect- 
ed would act most nearly like that of a wing with, full-snan 
flaps. Similarly, the tip section would behave as if it were 
on a wing with no flap deflection. Now the load distribution for 
the wing with partial span flaps can be calculated as long as 
each section is operating in the linear range of its lift curve. 

A typical distribution, denoted by ci g, is shown in Figure 5. 

A lift distribution is now calculated for the same wing 
with no flap deflection and the ordinates scaled so that the wing 
has a lift coefficient of unity. Such a distribution, denoted 
by c; ( is also shown in Figure 5. 

This distribution is multiplied by a suitable constant 
factor k 1 to give the same value at the root as the distribution 
with the flap deflected and by another factor k2 to give the 
same value as the flap distribution at the outermost station 
used in the computations. 

The differences between c l 8 and k, c* i inboard of the flap 
end and between c t 8 and kgc^ , outboard of the flap end, are 
divided by the difference between kj c/t and k2 c taken at the 
end of the flap. Thus, the resulting values or factors, F can 
be expressed as follows: 


c / 8~ k i o/i 
c/ | ( k| - k 2 ) 



(41) 


and for > J^_* 

b b 
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c/S-k 2 c / 1 

c /l (k| — k 2 ) 


(42) 


F =• 


These factors are relatively insensitive to lift-curve slope 
and may be used in the non-linear portions of the lift curves. 

In calculating the above mentioned distributions, any 
convenient values of the discontinuity S and wing angle of attack 
may be used. Wing twist is omitted in these computations. 

Having obtained the factors F the values of maximum lift coeffi- 
cient are then altered according to: 

max = ( c ^max)o + F(A c /max) (43) 


where ( c ?max)o is the two-dimensional value of section maximum 
lift coefficient and A max is the increment in c* m ax due 
to flap deflection at the end of the flap. 

Figure 6 illustrates the method of obtaining the final 
corrected lift curves. At each angle of attack (a 0 ) the two- 
dimensional uncorrected value of lift coefficient ( c/o) , as 
obtained from the airfoil data, is factored to give the corre- 
sponding value of lift coefficient corrected for proximity to 
the flap end thus, 


c £ cor = c /o 


c / max 
( c £ max)o 


(44) 


This value of lift coefficient would correspond to a value 
of angle of attack a 1 on the two-dimensional lift curve. Due 
to the fact that each wing section is not operating in a true 
two-dimensional flow, but forms part of a finite- span wing a 
correction must be made to account for the fact that the air 
can flow around the wing tips. This extra degree of freedom 
implies that a section on a finite- span wing must operate at a 
higher angle of attack than the same section in strictly two- 
dimensional flow in order to achieve the same value of lift 
coefficient. This correction was first derived by Jones, Ref- 
erence 48, in the form of the edge velocity factor given by: 


ae - g/o 
a 0 - a/ o 


(45) 


Jones originally derived the correction factor for a wing 
of elliptic planform and expressed it as. 
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Figure 6. Illustration of Method for Correcting Two-dimensional 

Section Data. 
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wing semi-perimeter 
h " wing span 


( 46 ) 


Reference 2 gives a more accurate expression, 
to wings of other than elliptical planform, in the 


applicable 

form, 




_ 4 _ 

/R2 


(47) 


Thus, the effective angle of attack, ae , corresponding to the 
two-dimensional angle of attack, a 1 , is given by, 


P - a e - Q; o 
a' - ai o 

and since a 1 is related to the original value do by 


o' -Q/o . c / mox 

“o-^o ( c * max)o 


(48) 


(49) 


then the correction which must be made to ct 0 to account for 
both edge-velocity and flap proximity is obtained as 


a e-a/ 0 _ c / mox 

a Q - Q/ 0 ( c / max)o 


(50) 


Figure 6 shows the final corrected lift curves as compared 
to the original two-dimensional data for two stations, one 
inboard and the other outboard of the flap end. 

The above method of correcting the two dimensional section 
data to account for flap discontinuity is substantiated 
experimentally by Reference 2. 

It should be noted that the Jones edge velocity correction 
applies regardless of whether the flap is deflected or not. 

3.2.5 Calculation of Overall Wing Characteristics 

Once the spanwise lift distribution of a wing has been 
calculated the determination of drag and pitching moment 
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coefficients is a simple matter since these quantities depend 
on the lift distribution. 

By the use of Simpson’s Rule, the integrated values of 
lift, induced drag, profile drag and pitching moment coefficients 
are given by: 


Cl ■$■)* 1 V- 

m=i 

c 0 ,=(i) 2 * i; 


(51) 


(52) 

(53) 


(R b 
c m ' b 2 c' 




m 


where 


c m = c m c / 4 ~ 



cos ( a B - a/ ) + c do sin ( a B - ) 


sin (a B -a/ ) - c do cos (a B -a;) 


(54) 


(55) 


In equation (55) the section lift and drag forces are assumed 
to act through the quarter-chord points. 

The multipliers used in equation (51) through (54) are 
defined by: 



Step-by-step computational procedures based on the preceding 
analysis are presented in the following section. 
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SECTION 4 
COMPUTER PROGRAM 

The solution of the mathematical model described in Section 
3 was accomplished by means of a specially developed digital 
computer program. This section presents a detailed description 
of the computational procedures, the required section data, and 
the implementation of the computer program. The results of 
sample calculations are also presented. 

4 . 1 COMPUTATIONAL PROCEDURES 


The following is a mathematical description of the individual 
steps used in the computation of the wing spanwise load distrib- 
utions. The sequence of the calculations is essentially the 
same as that set up in the computer program which is described 
later. The computational procedures cover three separate cases; 
wings with no flaps, wings with full- span flaps and wings with 
part-span deflected flaps. 

4.1.1 Computation of Basic Parameters 

(a) Calculate the following geometric quantities 


Yo=By,_ Hi 
A 2 

(57) 

e = V A 2 -B 2 

(58) 


4Yo (Vc)R( c R/b) 
7T A B 


(59) 


, 2y 0 2B 1 2 h 2 A 1 2 e 1 

where Yo.B.H.A.e are, respectively, — 

(b) Calculate a number of points y 1 on the exposed wing 
span for even increments of Q using 


Y 1 = Yo+ ( I -Yo ) cos 8 


(60) 
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where 


9 - cos 1 


= cos 


-i 


c-M-) 


(c) Compute the average (non-dimensional) distances 
the points Y 1 from the focii of the elliptic fuselage using 


= ~T [y Yi2 + (H-e) 2 + yyi2 + (H+e)2 ] 


where for the wing tip 


°T~ 2 I + (H-e) 2 + v^l~+(H + e) 2 j 


(d) Transform from the u -plane to the u -plane using 
following relationships for an elliptical fuselage. 


= aM' 


B Qt 


'n/q J - < 




and 


Y. 

b A B (. y a , 2 _ e 2 J 


If the fuselage has a circular cross-section compute 
b _ 

b 1 + H 2 
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(61) 

of 

(62) 

(63) 
the 

(64) 

(65) 

( 66 ) 


I 



and 


Y'= Y' 


A 2 

Y'2 + H 2 


(67) 


where A is the cross-sectional radius (non-dimensional), 

(e) Calculate a new set of points Y on the transformed 
span for even increments in the spanwise variable q as follows 


Y = cos Q (68) 

where 



(f) Using the relationship of Y 1 to Y 1 obtained above, 
interpolate to find the points on the physical span, Y , 
corresponding to Y. 

(g) For an elliptical fuselage cross-section obtain values 

for 


R 


du 

du 


M*- 


VZT 


I +- 


e 2 Y 2 
( a 2 _ e 2 }2 


(69) 


where 


a = \ (^/ Y 2 +(H-e ) 2 + 7Y 2 +(Hte) 2 ^ 


(70) 


or, if the cross-section is circular, use 


R du A 2 (Y 2 -H 2 ) 

- du * ( Y 2 + H 2 ) 2 


(71) 


44 



(h) If the wing has a deflected part-span flap calculate 
the location of the end of the flap in the u -plane using 


v*. b Y* 
Y ' "E (A-B) 



Bo* | 

v /a*2-e2 J 


( 72 ) 


where 


a * = z <\/y * 2 + (H e ) 2 + VY* 2 +(H + e ) 2 


( 73 ) 


(i) For a wing with linear taper in both chord and 
thickness from fuselage side to wing tip calculate 




'lYI-Yo 1 

. I -Yo J 


( 74 ) 


(t/c) . / x (t/c)T \r lYI-Yo " 

(t/c) D V A (t/c) R-'L 1 ~Yo_ 

K 


K 


K r + ( K t K r ) 


IYI - Yo ' 
I -Yo. 


( 75 ) 


( 76 ) 


g-r\ IYI - Y 

c/Cr I - Y 


( 77 ) 


and 


R e '( °/ c R ) 
e ~ c V c R 


( 78 ) 
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where Re 1 is the flight Reynolds number based on the exposed 
wing mean aerodynamic chord c 1 which is given by 

c 1 _ 2 I + X + X 2 1 (79) 

C R 3 [ I + X J 

Note that the geometric wing twist, e , is non-linear and that 
the spanwise distribution of wing section camber level is taken 
to be linear. 

For a wing planform which, is not trapezoidal the foregoing 
quantities may be determined from a drawing or from a special 
calculation. 

(j) Finally, calculate the multipliers /3mk , rj m , and the 
matrix of coefficients Kij using equations (23) (56) and (30) 
respectively. Also, for part-span deflected flaps obtain c/gc/BS 
and ack /8 from equations (38 ) and (39). 

4.1.2 Wings With No Flaps or With Full-Span Deflected Flaps 

If the wing has deflected or undeflected full-span flaps 
the computations are simplified since no corrections are required 
for spanwise discontinuities in induced angle of attack. In 
this case the computational procedure is as follows: 

(a) For the desired body angle of attack calculate the 
distribution of geometric angle of attack using 

a g = a R + e + a B T [o-^-l]+a B (80) 


(b) Obtain the corresponding c/ and a 0 values by interpola- 
tion in the two-dimensional section data at the proper Reynolds 
number, thickness-chord ratio and camber level. 

(c) Calculate an initial approximation to the distribution 
of the loading on the transformed wing using: 


c/c / /R 
b = C ( V /R + 1'8 


cX / 'cr\ ( b h 

A Cr A bA b/ 




(81) 


(d) Calculate the corresponding values of induced angle 
of attack <2/ using equation (21) and determine the effective wing 
angles of attack in the real plane from equation (11). 
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T 


(e) Now compute the equivalent angles of attack, ao , for 
use with the two-dimensional section data thus 


a o 




a e -a /o ( I — E) 



E 


(82) 


where E is given by equation (47). 

(f) Using the section data obtain the values of lift 
coefficient corresponding to ao and then calculate new values 
of c / c/b . 

(g) Compare these values to the guessed values of c/b 

and if agreement is not sufficiently close obtain a set of values 
of Aj by applying the matrix K t j to the differences between 
the calculated and guessed values as shown in equation^ ( 30b ) . 

Add the results to the original guessed values of cf c/b and 
obtain a new and better set of approximate values. Repeat steps 
(d) through (g) until satisfactory agreement is reached. 

(h.) Having determined the lift distribution obtain section 
values of profile-drag coefficient, induced-drag coefficient 
and pitching moment coefficient. 

(i) Finally, using equations (51) through ( 56 ) 9 calculate 
the values of Ci_ , Co 0 , CDj , and Cm . 

4.1.3 Wings With Part-Span Deflected Flaps 

(a) Determine the lift distribution for the plain untwisted, 
unflapped wing-body combination at an angle of attack which is 
within the linear range of the section lift curves. Scale the 
resulting lift distribution to give a value of wing lift coeffic- 
ient equal to unity. Denote this distribution by c/j . 

(b) For the untwisted, flapped wing body combination, select 
a body angle of attack within the linear range of both the 
flapped and unflapped section lift-curves. Obtain the value of 
section lift coefficient, c/ R , for the root section from the 
airfoil data. 

(c) Calculate the initial approximation to the lift distrib- 
ution in the transformed plane using 


c i c 
b 




< Y < Y ' 


(83) 
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I 


and 


etc j_ 
5 " 2 


c ‘" GSr) © 


(+){'V^W}* 


'£Y< I 


(84) 


(d) Taking the value of 8 as the difference between the 
flapped and unflapped zero lift angles of the section at the 
end of the flap, determine the resultant angles of attack, cte , at 
the real wing using equations (11), (36), (37), (38), (39). Correct 
these angles for edge velocity and then use the section data to 
obtain new values for the lift distribution. Repeat this itera- 
tive cycle until a convergence is obtained in the lift distrib- 
ution. Denote this distribution by §. 

(e) Calculate the correction factors F from equations (41) 
and ( 42 ) . 

(f) Using equation (81) obtain the first approximation 

of the final lift distribution and compute the value of section 
lift coefficient at the flap end, c/*, thus 


c /* 


(-¥-) <rB (4-) (-£) 


(85) 


(g) Determine the uncorrected value of lift coefficient at 
the end of the flap using 


c ' u 


c;* 

F F 


( 86 ) 


where FF_is the value of I + F A c /* max / ( c f max)o taken at the flap 
side of 7 * . 

(h) Obtain the angle of attack value, a 0 , which corresponds 
to c/ u , from the flapped section data and correct it to obtain 
the equivalent angle of attack, ae g thus 


a e § = E xFF (<*o - a/o) + a /o 


(87) 


(i) In the same way find the equivalent angle of attack, 
a /8 = o for the unflapped side of 7*. 



n 


f (j) The first approximate value of the discontinuity in 

angle of attack 8 is then given by 

S = a e g = 0 - aeg (88) 

(k) Again, use equations (11), (36), (37), (38), and (39) 
together with equation (88) to obtain values of the resultant 
angles of attack, a e , at each station on the real wing. 

(l) The corresponding uncorrected angles of attack are 
then found using 


Qo 


a e - Q/ o ( I - E ( c / max)o \ 
f c / max 
( c / max)o 


(89) 


(m) Use the two-dimensional section lift data to obtain 
values of c/ 0 at these angles of attack and factor them to give 
the correct values of lift coefficient, thus 


c ( z c /o 


H max 
( c / max)o 


(90) 


(n) Compare this calculated lift distribution to the 
approximate distribution and repeat the iterative process until 
the convergence is achieved within the specified limits. 

(o) Finally, obtain the overall wing-body characteristics 
using equation (51) through (56). 

4.2 SECTION DATA CHARACTERISTICS 


The accuracy of wing spanwise load distributions computed 
by this program largely depends on the quantity and quality of 
the available section data characteristics. The section charac- 
teristics required in the computations are: The two-dimensional 

section lift-coefficient versus section angle of attack, the 
section profile-drag and the quarter-chord pitching-moment 
coefficients versus section lift coefficient. Ideally, the data 
should be available for as broad a range of Reynolds numbers 
and thickness-chord ratios as is required to cover the range of 
these parameters expected in actual flight. Otherwise, linear 
interpolation and extrapolation of the existing section data are 
required to perform the computations at the values of Reynolds 
number and/or thickness-chord ratio outside the available range. 
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The most reliable airfoil section data which are now 
available are summarized in References 5 and 6. Reference 5 
presents the data obtained for Reynolds numbers between 3x10° 
to 9x10° while Reference 6 extends these measurements to values 
of Reynolds numbers down to 0.7x10^. 

4.2.1 Correction of the Section Data 

The senior author of Reference 6 has indicated that errors 
exist in the values of lift curve slope for all airfoils tested 
at Reynolds numbers of less than 3 . 0 x 106 . These errors are 
attributed to an angle of attack change associated with air leak- 
age at the intersection of the test airfoil and the wind tunnel 
walls. The evidence of the inconsistencies at low Reynolds 
numbers can be noted by examining the behavior of the values 
of section lift-curve slope as a function of Reynolds number, 
shown in Figure 17 of Reference 6. Since the values of lift 
coefficients obtained from force measurements are known to be 
correct the irregularities in the trends of c/ a versus Reynolds 
number can only be associated with the errors in airfoil angle 
of attack. 

Therefore, a correction procedure based on extrapolation 
of the high Reynolds number lift-curve slope data was employed 
to generate approximate section data for low Reynolds numbers. 
This procedure which was approved by the senior author of Refer- 
ence 6 is described below. 

For each airfoil section the variation of lift-curve slope 
with Reynolds number was plotted as shown in Figure 7. A 
straight-line extrapolation of the best fit to the data for 
Reynolds numbers of 3xlC)6 and greater was used to obtain values 
of lift-curve slope at low Reynolds numbers. These values were 
then used to correct the quoted angles of attack at constant 
values of lift coefficient. In making the extrapolations, 
account was taken of the variation of lift-curve slope with 
thickness-to-chord ratio at constant Reynolds number as shown 
in Figure 8. Furthermore, an effort was made to ensure that the 
extrapolation trends were compatible with the trends of the 
high Reynolds number data. The magnitude of the correction for 
a typical section is illustrated in Figure 9. Table I indicates 
which airfoil section data required this correction. 

4.2.2 Preparation of the Section Data for Computer 

The airfoil section data, corrected as required using the 
previously discussed procedures, are prepared in the form of 
special tables, which are used in the computer. For a given 
airfoil these data represent the relationships of c/ vs. a , c<j 0 
vs. c / , and c m c/4 vs. C/ for constant values of Reynolds number. 
Typical tabulation of these section characteristics suitable 
for the use in the computer are shown in Figure 10. 
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©Data from NACA TN 1945 


NAGA 64409 


NACA 64412 


NAGA 64415 


Figure 7. Extrapolations of Lift Curve Slopes at Low Reynolds Number 
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Figure 8. Variation of Section Lift- Curve Slope with 

Thickness-Chord Ratio at Constant Reynolds number 
NACA 644 Sections 
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Figure 9 - Corrected Lift Curves for NACA 64-421 Airfoil at Low Reynolds Numbers 















n L± 


Table of c i Values 



-90. (f 0 


1.5 

2.0 

3.0 

0 

-0.40 

0 

-0.50 

0 

-0.60 

l! 42 

l.*48 

U47 

l.*28 

10.0 

l.*30 

10.0 

l.*26 

10.0 

1.44 

1.49 

1.51 

13.4 

13.4 

13.4 




Table of Values 
Reynolds Number x!0“^ 


2.0 

3.0 

2.000 

2.000 

0.009 

0.007 

0.007 

0.007 

0.020 

0.019 

2.000 

2.000 


Table of c m Values 

\ Cf Reynolds Number xlO"^ 

/ ' 0'.'7 |T.0l' lT | 2.0 P‘370 

I ' pT -10.0* 0. 0. 0. 0. 0. 

Re 3 -0.68 -0.09 -0.09 -0.09 -0.09 -0.09 

Re 2 0.*99 -0i09 -o!o 8 -o! 09 -o!o9 -o!o9 

' / Re n ■••••• 

L 1.64 -0.06 -0.06 -0.06 -0.06 -0.06 

10*0* 0. 0. 0. 0. 0^ 

0 . 0. 0. 0. 0 . 

0 . | 0. [ 0. 0. I 0. 

* For computer continuity only. 
Figure 10. - Method of Tabulation of Section Characteristics. 



The tabulation of section lift data is considered to be 
most critical as compared to the two sets of data cited above. 

For a given airfoil thickness-chord ratio, values of lift 
coefficient are carefully read at selected angles of attack for 
each Reynolds number so as to best define the lift-curves 
especially in the vicinity of C; max.. Large scale plots of 
c/ vs . a (e.g. NACA originals) are recommended for this purpose 
for improved accuracy. 

The values of lift coefficient thus obtained are entered 
in an array as shown in Figure 10 with the first and last 
entries of angle of attack being -90° and +90° corresponding 
to lift coefficient values of zero and 10.0 respectively. The 
latter value of section lift coefficient is entered only for 
the purpose of the computer in order to avoid computational 
discontinuity in the case of wings with part- span fully deflected 
flaps. The bottom two rows of the array contain the values of 
maximum lift coefficient ( c /max) and the corresponding angles 
of attack at which these first occur ( a max). This procedure 
is repeated for each value of airfoil thickness-chord ratio for 
which data is available. 

Similarly, using the plots of c do vs. c l and c mc/4 vs. c / 
the drag and pitching moment tables are prepared and tabulated 
as shown in Figure 10. In these arrays of data the last two 
rows are entered as zeros. 

Thus, a library of section characteristics tables is prepar- 
ed. The data is key-punched and stored ready for use in the 
computer. The families of airfoils for which data has been 
prepared are summarized in Table I. 

4.2.3 Airfoil Data Table Look-up Procedures 


The basic airfoil data are read into the computer row-wise 
and stored on tape. During execution of the program the data 
is read from tape into core in the form of files. The arrange- 
ment and designation of the files is shown diagrammatically 
in Figure 11. Each file must contain at least two, but not more 
than five, levels or values of thickness-chord ratio. Each level, 
or each table of data, is limited to 25 rows and 12 columns. 

This file size is adequate for most purposes, but may be 
increased if desired by changing the specifications relating 
to file size in the program. 

Linear interpolation is used throughout and is performed 
first for the required value of Reynolds number then thickness- 
chord ratio and finally for the given value of camber. When the 
required angle of attack is outside the range of the maximum 
angles of attack listed in the table for the two bracketing 
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c c vs. a 




: , d vs. c / 


Subrout me BRI DG 
interpolates 
between 
camber levels 


Subroutine 'ARC' 
interpolates 
between 

thickness levels 



(t/c)]_ 
C t/c ) 2 
(t/c) 3 
(t/c) 4 
( t/c ) 5 


File 5 cc vs. a 



File 10 c d vs. c/ 



vs. c i 


C 


Subroutine 'LOOK' 
interpolates 
within 
each level 


File 4 c/ vs. a 
with flap deflected 


File 15 c m vs. 



File 20 vs. a 
with flap deflected 



Airfoil data for root Airfoil data for tip 

camber level eg. 230XX camber level eg. 430XX 

Figure 11. - Schematic Representation of Section Data Storage 
in the Computer. 




values of the parameters under consideration, simple linear 
interpolation formulae are utilized with no special order of 
interpolation, A special computational procedure is required 
to obtain values of lift coefficient when the angle of attack 
falls between the tabulated values of C*max. This procedure 
is briefly described in the following paragraphs as an aid 
to understanding the program listing. 

Using Figures 12(a) and 12(b), the values of a max and 
c * max at the required value of Reynolds number, Re , are determined 
for each of the bounding values of thickness-chord ratio, (t/c)i 
and(t/c )2 . This yields values of lift coefficients at the 
points nn| ? and nri 2 , respectively. Depending on whether the given 
value of angle of attack, cl y is greater or less than the values 
of a mox, appropriate values of c/ (denoted by x i , *2 ) are obtained 
for each bracketing value of thickness-chord ratio. 

Knowing the values at points mj andm 2 , new values of c /max 
and a mox are calculated corresponding to the required value of 
thickness-chord ratio and are denoted by the point . Using 
the previously computed values at points X| (orX 2 ) together 
with the values at m 3 , the above interpolation cycle is repeated 
to obtain the value of at the required thickness ratio. 

This value of C/ now corresponds to the required values of 
Reynolds number and thickness-chord ratio, but applies only to 
the camber level associated with the tip airfoil series. The 
aerodynamic data for the root airfoil series is then called into 
core and the process is repeated. The final interpolation shown 
in Figure 12(d) is performed for the correct value of camber 
level in the same way as that used for thickness-chord ratio. 

4.3 DESCRIPTION OF THE COMPUTER PROGRAM 


The computational procedures described in Section 4.1 have 
been programmed for use on the NASA CDC 6600 series computer 
located at Langley Research Center. The program is written in 
the Fortran IV machine language. An internal listing of the 
program is presented in Appendix A and a block diagram illustra- 
ting the major logic features is presented in Figure 13. 

The program is initiated by reading in the basic configura- 
tion parameters punched on cards as indicated in Figure 14. 

The values of either aerodynamic or geometric twist may be input. 
Whichever twist is specified the column reserved for the other 
twist must contain 100. Columns 1 to 11 on card #2 contain the 
maximum allowable difference between the calculated and 
approximate values of c / C// b at the end of each iteration. A 
value of .001 was used in the computations presented later 
in this report. 
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FIGURE 13 COMPUTER PROGRAM BLOCK DIAGRAM 
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Figure 14. - Schematic Representation of the Computer Input Cards. 
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Normally the distribution of such planform parameters as 
chord, twist, local Reynolds number, etc. are calculated but 
provision is made to read these in if the planform is not 
trapezoidal. The program performs the computations using 10 
control points per wing semi-span. This number of stations is 
usually adequate for most purposes, however, if it is desired to 
increase (or decrease) the number of points the program can 
easily be modified. 

Having calculated, or read, the geometric parameters and 
computed the parameters governing fuselage transformation (if a 
fuselage is present), the values of the multipliers /3 mk and r){ 
are computed and stored. If the calculation is for a wing with a 
deflected part-span flap the values of the parameters associated 
with a spanwise discontinuity in angle of attack are computed. 

The coefficients of the matrix Kfj , given by equation (30b), are 
now computed using a matrix inversion subroutine to obtain success- 
ive approximations to the lift distribution. If the calculation 
is for a flapped wing the program branches to compute the two 
load distributions required to obtain the factors used in 
altering the two-dimensional section lift data. Having obtained 
these factors the basic program iterative loop is entered and 
executed until a lift distribution is calculated which agrees 
with the guessed distribution to within the required tolerance. 

The computed distributions of section lift coefficient, etc., 
and the corresponding integrated values are immediately printed 
out for each case run. 

If, at any point on the wing, the stall is detected, i.e. the 
computed value of effective section angle of attack exceeds that 
for maximum lift, the computations are then repeated for a value 
of angle of attack half-way between the last two values. If 
at this intermediate value the wing unstalls the program increases 
the angle of attack by 0.2 degrees until stall again occurs or 
until the increments total 1°. Similarly, if stall is still 
detected at the intermediate value, the angle of attack is 
decreased until the wing is unstalled or until the total decrease 
equals 1°. 

Alternatively, if, for a given body angle of attack wing stall 
is not obtained, the computations are performed for the next 
inputed value of body incidence until the value of a b =99 is 
encountered whereupon the calculations are stopped. A typical 
print out of the final results is presented in Table II. 

The computer program described herein can be used to predict 
the distributions and the integrated values of lift, drag and 
pitching moment coefficients for wings of trapezoidal planform 
with zero sweep in the presence of a fuselage at all angles of 
attack up to and including stall. In addition, these computations 
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642 SERIES FLAP CASE 


SPANWISE STATIONS 


( l) 

5*57711 503E-0 1 

( 2) 

9 * 6 1459 121 E-0 1 

( 31 

9. 02092175 E- 01 

( 4) 

8.2118 5699E-01 

( 5) 

7. 2 0940345 E -01 

( 6) 

6 • C4 2 191 1 5E- 0 l 

( 7) 

4.7471 2490E -0 1 

( 8) 

3*37573147 E— 01 

( 9) 

2 • 02203541E-01 

(10) 

1.00000444E-01 

( II ) 

-2.02203541E-01 

( 12) 

-3.37573147E-01 

(13) 

-4.74712 490E -01 

( 14) 

-6.04219115E-01 

(15) 

-7. 20940 345 E-0 1 

( 1 6 ) 

-6*21 1 85699E — 01 

(17) 

-9.02092 175E-01 

(18) 

—9. 61459121 E- 01 

( 19) 

-9 .9771 1 503E— 01 

( 



ALPHA MAX 









< I) 

1 • 6088983 1 E+0 1 

( 2) 

1 • 6 16 3520 7E + 0 1 

( 3) 

1 • 6322 7454E + 0 1 

( 4) 

1*6641 726 3E +0 1 

( 5) 

1.7260341 IE *-01 

I 6) 

1* 8 5084020E *0 l 

( 7) 

2.37034870E+01 

( 8) 

1 .25746 196E+0 1 

( 9) 

1 • 39 74535 OE +01 

(10) 

1 .43900 028E *01 

(III 

1 • 39745350E + 0 l 

(12) 

1 . 2 5746 196E + 0 1 

(131 

2.37034870E+01 

( 14) 

1 . 85C 84020E+01 

(15) 

1.7260341 1E+01 

(16) 

1.66417263E+01 

(17) 

1.6322 7454E *-01 

( 18) 

1 • 61 635207E*0 1 

(19) 

1 • 6088983 IE **01 

( 



CL MAX 









( 1) 

1.55000003EfQ0 

( 2) 

1 • 556 5 0944E *00 

( 3) 

1.57042592E+00 

( 4) 

1 • 5982992 9E +00 

( 5) 

1.65222 564E+00 

( 6) 

1 • 76 061 569E +00 

( 7) 

2.12708990E+00 

( 8) 

2* 5583496 1E + 00 

( 9) 

2*6558441 OE+OO 

(10) 

•2.73499988E+00 

cm 

2 • 6958441 OE+OO 

1 12) 

2*55834 96 1E*00 

( 13) 

2 • 12708990E *00 

(14) 

l . 760 6 1 56 9E *-00 

(15) 

1 .65222 564E *00 

(16) 

1 • 5982992 9E + 00 

(17) 

1.57042592E+00 

( 18) 

1.5565C944E+00 

( 19) 

1 • 5500C003E+00 

( 



THICKNESS / 

CHORD 

DISTRIBUTION 







( 1) 

1*20152 566E-0 1 

( 2) 

1.22569392E-01 

( 3) 

1 • 2652 71 88E-0 1 

( 4) 

1. 3192C953E-01 

( 5) 

1. 38603977 E- 01 

( 6) 

1 • 46 385392E— 0 1 

( 7) 

1*5501 916 7E— 0 1 

( 8) 

1*6416 1 790E— 0 1 

( 9) 

1*7318643 IE-01 

(10) 

1 . 79999970E-01 

(11) 

1*73 18643 IE— 0 1 

(12) 

1*64161790 E— 01 

( 13) 

1 • 550 1 9 L67E-0 1 

( 14) 

1.46385392E— 01 

(15) 

1.38603977E-01 

(16) 

1*31 920953E- 0 1 

(17) 

1 . 2652 7 1 88 E-0 1 

(18) 

1 *22569392 E-0 1 

( 19) 

1.201 5 2566E-01 

( 



SECTION REYNCLDS 

NUMBER 







( 1) 

6. 00060006E*00 

( 2) 

6.00060006E+00 

( 3) 

6.00060006E+00 

( 4) 

6.0C0600C6E *00 

( 5) 

6.00060006E*00 

( 6) 

6. 00060006E+00 

( 7) 

6.000600C6E+00 

( 3) 

6.00060006E+00 

( 9) 

6.000600C6E+00 

(10) 

6.00060006E+00 

(11) 

6.0C060006E+00 

(12) 

6 • 0006 0 006 E + 00 

(13) 

6.00060006E+00 

( 14) 

6.0006C006E+00 

(15) 

6. 0 006 00 06 E *00 

(16) 

fc .00060006E *-00 

(17) 

6*00060006E*00 

(18) 

6.0 00 600 06 E *-00 

( 19) 

6 . 0006000 6E +00 

( 



CHORD DISTRIBUTION 







( 1) 

1.00000000E**00 

( 2) 

l.OOOOOOOOE+OO 

( 3) 

l.OOOOOOOOE+OO 

( 4) 

l.COOOOOOOE+OO 

( 5) 

l.OOOOOOOOE+OO 

( 6) 

1.00000000E+00 

( 7) 

l.OOOOOOOOE+OO 

( 8) 

1.000000006+00 

( 9) 

l.OOOOGOOOE+OO 

(10) 

1.00000000E*00 

(11) 

l.OOOOOOOOEfOO 

(12) 

l.OOOOOOOOEfOO 

( 13) 

1.00000000E*-00 

(14) 

l.CCOOOOOOE+OO 

(15) 

1.00000000E*00 

(16) 

l.OOOOOOOOE+OO 

(17) 

l.OOOOOOOOE+OO 

( 13) 

l.OOOOOOOOE+OO 

( 19) 

l.OCOCCOCCE+OO 

( 



GEOMETRIC TWIST 








( 1) 

-4.92881352E+00 

( 2) 

-4.72977270E+00 

( 3) 

-4.40382321E+00 

( 4) 

-3.95561265E +00 

( 5) 

-3.40922351E+00 

( 6) 

-2*76 837490E+00 

( 7) 

-2.05732909E*00 

( 8) 

-l*30437645E*-00 

( 9) 

-5*611 4040 5E -01 

( 10) 

-2.43856 132 E- 06 

(11) 

— 5* 61 140405E-0 L 

(12) 

-1.30437645E+00 

( 13) 

-2.05732909E+00 

( 14) 

-2.76837490E+00 

(15) 

-3.40922351E+00 

( 16) 

-3.55961265E+00 

(17) 

-4*40382321E+00 

(18) 

-4.72977270E+00 

( 19) 

-4. 9288135 2E+ 00 

( 
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642 SERIES FLAP CASE 


./, 


/.i/. 


BODY ANGLE OF ATTACK, DEG. 

BODY HEIGHT /SPAN 

ASPECT RATIO 

WING BODY INCIDENCE, DEG . 
ROOT THICKNESS CHORD . . . 
NUMBER OF SPANWISE STATIONS 
FLAP SPAN / WING SPAN. . . 
FLAP SETTING, DEG 


16.00 

.10 

6.00 

0.00 

.18 

20.00 

.45 

60.00 


COORDINATES OF MOMENT REFERENCE POINT 


VALUE OF DISCRIMINANT. . 
BODY WIDTH / SPAN. . . . 
WING HEIGHT /SPAN. . . . 
TIP THICKNESS CHORD. . . 
GEOMETRIC TWIST, DEG . . 
AERODYNAMIC TWIST, DEG . 

TAPER RATIO 

REYNOLDS NUMBER 

X= 0.00 t- 


.001000 

.10 

0.00 

.12 

-5.00 

1.00 

6.00 


0.00 


f. ./../. 







DISTRIBUTION 

GF SECTION LIFT COEFFICIENT 






( n 

3.42250509E-01 

( 2) 

6.30955 94 OE— 01 

( 3) 

8. 65952551 E- 01 

( 4) 

1. 061 70561E +00 

( 5) 

1.2411 46 08E+00 

( 6) 

1.43723 53 7E +00 

( 7) 

1 • 8298 1 852E +00 

C 8) 

2. 23042330E+00 

( 9) 

2.4226873 7E +00 

(10) 

2. 543 5842 2E + 00 

( m 

2.42268686E+00 

(12) 

2 • 23042 07 4E + 00 

( 13) 

1. 829 85557E + 00 

( 14) 

1.43723226E+00 

(15) 

1.241 14578E+00 

( 16 1 

1 . C 6 1 70626E + 00 

( 17) 

8.65954479E-01 

(18) 

6. 30959 08 5 E- 01 

( 19) 

3.4225474 4E— 01 

( 



STALL MARGIN 

DISTRIBUTION 







c i) 

1 * 2 C 77495 26 + 00 

( 2) 

9. 255 53 49 8E— 01 

( 3) 

7. 044 73366 E-0 1 

( 4) 

5. 365936 79E-01 

( 5) 

4. 110795578-01 

( 6 1 

3.23 380 322E-0 1 

< 7) 

2. 972 7 1 378E-01 

( 8) 

3. 2792631 OE-O t 

( 9) 

2.7315 6732E-01 

(10) 

1.91415661 E— 01 

(11) 

2.73 15 7242 E-0 1 

< 12) 

3.2792 8 87 3E-0 1 

(13) 

2. 9723 43 34 E- 01 

( 14) 

3.2338 34 36E-01 

( 15) 

4. 110798598-01 

(16) 

5. 36 59302 3E-01 

< 17) 

7.G4471438E-01 

( 18) 

9.25550353E-01 

( 19) 

1 . 207 745 2 9E +00 

( 



BODY ANGLE OF 

ATTACK, DEG. . = 


16.00 LIFT COEFFICIENT . . . , 


= 1.666010 


TABLE II 


Concluded 



incorporate the effects of linear spanwise variation in section 
camber for unflapped wings. The effect of a deflected part-span 
flap on the lift distributions is included. However, drag and 
pitching moment results are not included for these cases because 
of the lack of section drag and pitching moment characteristics 
for flapped airfoils. 

The computer program described above constitutes a part 
of this report and is therefore available for public use through 
Computer Software Management and Information Center (COSMIC), 
as stated in Appendix B. 

4.4 RESULTS OF SAMPLE CALCULATIONS 


The computer program was used to calculate the character- 
istics of three unflapped wings and one wing having a 60% span 
split flap. These computer results were correlated with the 
corresponding theoretical predictions of References 1 and 2 and 
the available experimental data. Figures 15, 16, and 17 present 
computed lift, drag and pitching moment characteristics for 
unflapped wings of aspect ratio 8.04, 10.05 and 12.06, and Figure 
18, shows the variation of lift coefficient with angle of attack 
for a flapped wing of aspect ratio 9.02, having a 60% span, 

20% chord split flap. Included in these figures are typical 
spanwise distributions of section lift, drag and pitching moment 
coefficients. 

The correlation of the computer results versus theoretical 
characteristics of References 1 and 2 and the available exper- 
imental data are performed, where possible on the basis of the 
integrated values as well as the spanwise distributions of wing 
characteristics. Unfortunately, measurements of wing character- 
istics in the presence of a fuselage, especially of the span 
loading, are extremely difficult to obtain for the class of wing- 
fuselage combinations considered herein. As a result, there is 
a lack of suitable experimental data of the type required for 
the present comparisons. Some experimental data exists for swept 
wing body combinations but such data is not applicable for 
correlations with the straight wing results. 

Examining Figures 15 through 18, it can be noted that the 
computer results are in good agreement with the available exper- 
imental data. The slight differences between the present compu- 
tations and those of References 1 and 2 are attributed to hand 
versus machine fairing of the section characteristics. 
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Lift coefficient, 



Geometric Washout = 4.5° Root Section 4416 

Aspect Ratio =8.04 , Tip Section 4412 

Reynolds Number = 4.32x10° Taper Ratio = .4 


Figure 15. Experimental and Calculated Characteristics for a Wing 
of Aspect Ratio 8.04. 
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Figure 15. Concluded 
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Lift coefficient 



Drag coefficient, C n Angle of attack, „ Pitching-moment 

u B coefficient, C m 

Geometric Washout = 3.5° Root Section 4420 

Aspect Ratio = 10.05 , Tip Section 4412 

Reynolds Number = 3.49x10° Taper Ratio = .4 


Figure 16. Experimental and Calculated Characteristics for a Wing of 
Aspect Ratio 10.05 







Lift coefficient, Cl 



Drag coefficient, C D Angle of attack, «c B coeffiS^*^ 


Geometric Washout = 3,0° 
Aspect Ratio =12.06 
Reynolds Number = 2.87x10° 


Root Section 4424 
Tip Section 4412 
Taper Ratio = .4 


Figure 17. - Experimental and Calculated Characteristics for a Wing 
of Aspect Ratio 12.06. 



mc/4 Drag coefficient 




Figure 17 . Concluded 
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Ca) Calculated span load distribution for 10° angle of attack. 



(b) Comparison of experimental and calculated lift curve. 


Figure 18. Experimental and CalculatedICharacteristics for Wing 
with 60% Flap; Aspect Ratio 9.02; Taper Ratio 
0.4; Washout 2°. 
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SECTION 5 


PARAMETRIC INVESTIGATION 


The computer program was used to determine the effect of 
wing geometry on the spanwise lift distribution, the location of 
the initial stalling point and the corresponding value of the 
effective maximum wing lift coefficient, Cl max. A similar 
investigation based on a linearized formulation of the lifting- 
line theory is reported in Reference 49. The present results 
may be considered as an extension of this work in that a broader 
range of parameters was investigated and the more accurate non- 
linear theory was employed. 

5.1 RANGE OF PARAMETERS 


The range of parameters selected for this study is representa- 
tive of that applicable to present-day light aircraft. 

Computations were performed for wings utilizing the three basic 
airfoil series commonly encountered, NACA 64 series, 44 series 
and 230 series. These calculations were selected to show 
the effect of the three major geometric parameters, aspect ratio, 
taper ratio and section thickness, together with the influence 
of washout and linear camber variation from root to tip. Varia- 
tions in flight Reynolds number, based on the wing mean aerodynamic 
chord, were also investigated as were the influences of a part- 
span deflected flap and the presence of a fuselage. 

The variation of wing thickness and chord length was linear 
in all cases but the distribution of washout was non-linear in 
contrast to that used in Reference 49. The non-linear washout 
distribution chosen ensures straight-line leading and trailing 
edges. Washout was aerodynamic which is defined as the angle 
between the zero-lift lines of the root and tip sections. All 
wings were of trapezoidal planform without rounded tips, and zero 
sweep . 

The investigation encompassed 331 different configurations 
which are summarized in Table III. The configuration defined by an 
aspect ratio of 6, a Reynolds number of 6 x 106 an d a root and tip 
thickness chord ratio of 0.18 and 0.12 respectively was selected 
as a standard case for systematic variation of the parameters. 

These variations included the values of aspect ratios of 6, 8, 
and 10, aerodynamic washout of 0°, 2%°, 5°, and 7%°, taper ratios 
of 0.5, 0.75 and 1.0; root thickness ratios of 0.21, 0.18, 0.15, 
and 0.12; and tip thickness ratios of 0.12 and 0.15. 

The complete range of calculations were not performed for 
untapered wings having 7%° washout since this amount of twist 
would not normally be used on wings of this planform. However, 
to assist in establishing trends, a few computations were obtained 
for this washout. 
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TABLE III - SUMMARY OF CONFIGURATIONS STUDIED 


(a) Wings of NACA 44XX Airfoil Section 


Taper ] 
Ratio 

Aspect 

Ratio 


Root 

t/c 

i 

Tip t/c 

Aerodynamic 
Washout -Deg 

, • 

Re. No. 
X 10" 6 

Comments 

• 

.75 

.5! 

J 

6 8 

10 

.12 

.15 

.18 

. 21 ! 

! 

.12 .15 .18 

0 

2.5 

5.0 7 

.5 

3 6 9 

X 

X 

X i 

X 




X 

I 

X 

X 

X 

X 


X 


X 

X 

X j 

X 




X 

! 

X 

X 

X 

X 


X 

Aspect 

X 

X 

X 


X 



X 

i 

X 

X 

X 

X 


X 

Ratio 


X 

X 

X 


i 


X 

i 

X 




X 

X 

Variation 


X 

X 

X 




X 


X 




X 

X 



X 

X 


X 



X 


X 

1 



X 

X 


X 

X 

x i 

X 

i 




X I 

X 

X 

X 

X 

1 

X i 


X 

X 

X 

X 

1 


X 



X 

X 

X 

X 


X 


X 

X 

X 

X 


X 



1 

i 

x 1 

X 

X 

X 


X 

Root (t/c) 


X 

X 

X 


1 



1 

X , 

x 1 

1 



X 

X 

Variation 


X 

X 

I x 



X 


1 

X 




X 

X 



X 


x 


X 



i 

X 




X 

X 


X 

X 

X 

i X 




X 

i 

X 

X 

X 

X 


X 


X 

X 

X 

X 




X 


X 




X 

X 


X 

X 

X 

X 




X 


X 

X 




X 

Tip (t/c) 



X 

X 


i 


X 


X 

! 

X 



X 

Variation 

X 



X 


i 

1 


X 


X 




X 

X 




X 

X 


[ 


X 


X 




X 

X 


X 

X 

X 

X 




X 


X 

X 

X 

X 


X 


X 

X 

X 

X 




X 


X 

X 

X 

X 


X 

Reynolds N 

X 

X 

X 

X 




X 


X 




X 

X 

Variation 

X 

X 

X 

X 




X 


X 




X 

X 



Total number of configurations = 106 


TABLE III - SUMMARY OF CONFIGURATIONS STUDIED - Continued 
(b) Wings of NACA 230XX Airfoil Section 



Total number of configurations = 99 


TABLE III - SUMMARY OF CONFIGURATIONS STUDIED - Continued 


(c) Wings of NACA 642XX Airfoil Section 
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Total number of configurations = 99 



TABLE III - SUMMARY OF CONFIGURATIONS STUDIED - Concluded 

(d) Effects of Camber, Fuselage and Flaps on Wings of Various NACA Airfoil 
Sections. 
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Total number of configurations = 27 
Grand total of configurations studied = 331 


5 . 2 METHOD OF PRESENTATION OF RESULTS 


The effects of various design parameters on the stalling 
characteristics of unswept wing aircraft are herein presented 
in the form of plots showing the stall margin distributions, 
location and movement of stall boundaries and the maximum values 
of integrated wing lift coefficient. 

5.2.1 Stall Margin Distributions 

The stall margin, A C( , is defined as the difference between 
the maximum section lift coefficient and the section lift coeffi- 
cient when the stall first occurs on the wing. Figure 19 shows 
typical distributions of C/ ,C/max and A C / computed for the 
conditions indicated. The spanwise location of zero stall margin 
corresponds to the point of onset of stall and the rate of separa 
tion of the A C/ curve from the horizontal axis indicates the rate 
of stall propagation across the span. A wing is usually assumed 
to have sufficient stall margin if a value of AC/ =0.1 is indica- 
ted at the 70% semispan station. 

Since it would be impractical to present all three distrib- 
utions (i.e. C/ , Of max , and Ac/ ) for each of the 331 cases 
computed, in the manner presented in Figure 19, only the distrib- 
utions of stall margin are given. 

5.2.2 Stall Boundaries 


In general, the curve representing the spanwise lift distrib 
ution at the onset of wing stall is tangent to the section C/max 
curve at more than one point. This implies that a few adjacent 
wing sections (i.e. a portion of the wing span) can be stalled 
simultaneously at a given operating condition. The limiting 
values of spanwise locations encompassing these stalled sections 
are defined as the stall boundaries. Using stall margin distrib- 
utions, such as shown in Figure 19, the inner and outer boundar- 
ies of the stalled wing area can be estimated as the spanwise 
locations where AC/= 0.01. Furthermore, the movement of these 
stall boundaries affecting the growth and the propagation of 
wing stall areas can be expressed as a function of the basic 
design parameters. For example, using Figures 19(a) and 19(b) 
the stall areas for the unflapped and the flapped wing configura- 
tions considered are defined by the stall boundaries 50% to 68% 
and 10% to 40% of the wing span, respectively. 

5.2.3 Maximum Lift Coefficient 

Maximum lift coefficient is an important criterion in 
assessing wing stall characteristics as affected by the variation 
of basic design parameters. This lift coefficient is defined 
as the integrated value obtained from the spanwise lift distribu- 
tions when any one of the wing sections is stalled. Such stall 
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Section lift coefficient, c/ Section lift coefficient 



Fraction semispan 

(a) Unflapped wing - 230 series sections. 



(b) Flapped wing - 60% span flap - 642 series sections. 


Figure 19. - Typical Lift Distributions Along Wing Span 
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Stall margin, * c/ Stall margin, a c; 




is obtained when the local section lift coefficient 0/ equals 
the maximum value (C/max) of that section. This value depends 
primarily on section thickness-chord ratio, Reynolds number 
and the airfoil type. Typical variations of the maximum section 
lift coefficient (Qmax) with these parameters is shown in Figure 
20. The effects of various design parameters on integrated 
maximum wing lift coefficient are herein presented in the form 
of carpet plots. These plots facilitate linear readings of non- 
linear relationships between three variables and thus permit 
more accurate interpolation (or extrapolation) for intermediate 
values of the parameters. However, since this method of data 
presentation is not as commonly used as the familiar X-Y plot, an 
explanation of the use of carpet plots is presented in the 
following sub- sect ion. 

5.2.4 Carpet Plots 

The procedure for constructing and reading carpet plots can 
be best illustrated by the following example: The carpet plot 

shown in Figure 21 represents the variation of CLmax with Reynolds 
number and taper ratio for a wing employing 230 series sections 
and a constant value of washout of 7%°. Curve (abc) on this 
carpet is a conventional (x-y) plot of CLmax v Reynolds number 
for a constant value of taper ratio of 0.5. The indicated hor- 
izontal unit distance corresponds to a change in Reynolds number 
of 1x10^. Thus, point (a) corresponds to Reynolds number of 
3x10^ and point (b) which is 3 horizontal units from (a), corres- 
ponds to 6x10^ and point (c) to 9x10^. The curve (def) presents 
a similar plot for a taper ratio of 0.75. 

However, instead of plotting data for X = .75 on the same 
horizontal scale as that for curve (abc) a new scale is chosen 
whose origin is 2^ horizontal units to the right of the point 
(a), reflecting the change in taper ratio from 0.5 to 0.75. Thus, 
as regards changes in taper ratio, each horizontal unit corres- 
ponds to a change in X equal to 0.1. Similarly curve (ghi) for X 
= 1.0 is plotted 5 units to the right of point (a) or 2^ units 
to the right from point (b). After a new origin is selected 
for each taper ratio (i.e. point (a) for X= 0.5, point (d) 
for X = 0.75 and point (g) for X= 1.0), the plotting of the CL max 
vs. Re curves for each taper ratio is accomplished in the conven- 
tional way. If all points corresponding to the same value of CLmax 
Reynolds number are now joined, e.g. curves (adg) (beh) and (cfi), 
the resulting curves show the variation of Cl m ax with taper 
ratio for constant values of Reynolds number. 

Reading the carpet is as simple as constructing it. If 
the value of CLmax is required, say, for X = 0.62 and Re = 

5.3x10° it can be obtained as follows: For all three values 

of taper ratio locate and join points (j), (k) and (1) correspond- 

ing to Re = 5.3x106 . This is accomplished by moving 2.3 
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Figure 20. Continued - 44 series sections 




ma: 



Figure 21. - Variation of Cun ax with Reynolds 
Number and Taper Ratio. 
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horizontal units towards the right from point (a) to locate 
point ( j) on the X = 0.5 curve, 2.3 horizontal units from point 
(d) to locate point (k) on the X = 0.75 curve and the same 
number of units from point (g) to locate point (1) on the X = 1.0 
curve. Similarly, locate and join point (m) (n) and (o) moving 

1.2 horizontal units from points (a), (b) and (c) on the corres- 
ponding Reynolds number curves, respectively. The point of inter- 
section (p) of the curves (ihl) and (mno) yields the required 
value of CLmax = 1.445 for X = 0.62 and Re = 5.3x10^. 

5 . 3 COMPUTER RESULTS 


Presented in this section is a compilation of wing stall 
design charts obtained from the computer program. These charts 
can be used in preliminary design of unswept wing aircraft for 
determining optimum geometric parameters of a wing to yield good 
stalling characteristics. The effects of these geometric par- 
ameters on stalling behavior of straight wing aircraft are 
discussed in the following pages. 

5.3.1 Effect of Aspect Ratio 

Figures 22, 23, and 24 show the effect of aspect ratio on 
wing stall margins, stall boundaries and the values of maximum 
wing lift coefficient, respectively. Each of these figures is 
presented for three different families of airfoil sections, i. e. 
230, 44 and 642 series, and for a range of wing taper ratios 
between 0.5 and 1.0 and wing washout between 0° and 7%°. 

Examining Figure 22 it can be noted that an increase in 
aspect ratio results in a reduction of the stall margins at out- 
board wing stations. This effect is minimized at high values 
of washout and wing taper ratios. Furthermore, as can be seen 
from Figures 23 and 24, an increase in aspect ratio appears to 
have little effect on the spanwise location of stall boundaries 
and yields only a small increase in wing maximum lift coefficient. 

It can therefore be concluded that although wing aspect ratio 
is important from aircraft performance considerations, its 
effect on the stalling characteristics is very small. For this 
reason the remaining stall results which are presented in this 
section for aspect ratio of 6 only are considered to be typical 
and representative for the range of aspect ratios associated 
with present day light aircraft. 

5.3.2 Effect of Taper Ratio 

Figures 22 through 24 also show that taper ratio is one of 
the most dominant design parameters affecting wing stall charac- 
teristics. For any fixed aspect ratio and wing washout an 
increase in taper ratio from X = 0.5 to X = 1.0 (rectangular 
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Figure 24. Effect of Aspect Ratio and Taper Ratio on 0^^. 
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WASHOUT = 0 Deg. 


WASHOUT = 2% Deg. 


WASHOUT = 5 Deg. 


WASHOUT = 7h Deg, 
Root Section 4418 
Tip Section 4412 
Re = 6x10® 


Figure 24. Continued 
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WASHOUT = 0 Deg. 
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Root Section 64218 
Tip Section 64212 
Re =* 6x10° 


Figure 24. Concluded 
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wing) results in a substantial increase in stall margin on 
the outer portion of the wing. This increase in stall margin 
is associated with a shift of stall boundaries toward the inner 
portion of the wing. 

The advantage of the favorable stalling characteristics 
obtainable by increasing wing taper ratio is somewhat offset by 
a reduction of wing Cl max and hence, for a given wing loading, an 
increase in stalling speed. This reduction, which can be as 
high as 12% based on the Cl max values for a rectangular wing, 
depends on the airfoil sections and wing washout. 

5.3.3 Effect of Wing Washout 

Washout is often used to promote desirable stalling charac- 
teristics. For a given value of wing lift coefficient, twisting 
the tip section relative to the root section reduces the values 
of section lift coefficient in the tip regions and increases 
section lift coefficient at inboard stations. Since the distrib- 
ution of section maximum lift coefficient is unaffected by twist, 
the net result is to increase the stall margins in the tip region 
and to shift the stall boundaries inboard of the wing. These 
effects can be seen from Figures 22 and 23. 

Examining these figures it can be noted that for a taper 
ratio of 0.5 washout exerts relatively small influence on the 
stall margins at 70% semispan. As taper ratio increases, however, 
the effectiveness of washout becomes greater and is most effective 
for a rectangular wing (X = 1). However, a rectangular wing is 
unlikely to require washout since, as discussed previously, wings 
of this planform normally exhibit good stalling characteristics. 

The beneficial effect of washout appears to be offset by the 
adverse effect of the linearly diminishing Reynolds number. There- 
fore, the use of washout seems most justified for moderate values 
of wing taper ratio. 

Furthermore, Figure 24 indicates that for fixed values of 
aspect ratio and taper ratio, wing washout tends to reduce the 
value of CL max . The amount of this reduction depends on the wing 
sections and the value of washout used. 

In addition to the effects discussed above, washout may 
affect wing performance through changes in both wing profile 
and induced drag. The increase in induced drag coefficient 
due to washout above that for an untwisted wing is presented in 
Figure 25 for two values of aspect ratio 6 and 10. The data is 
for a wing having 230 series sections with root and tip thickness- 
chord ratios of .18 and .12, respectively. 

As can be seen, at low values of lift coefficient the use 
of washout causes an increase in induced drag while at high values 
of Cl washout reduces induced drag. Increased induced drag levels 
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Figure 25. - Increment of Induced-Drag Coefficient Due to Washout, 

230 series airfoil section - aspect ratio = 6 , Re = 6x10^ 
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at low lift coefficients may impair cruise performance depending 
on the amount of the increase relative to the overall drag 
coefficient of the airplane. Changes in profile drag due to 
washout were found to be sufficiently small that their effects 
on wing performance can be neglected. 

5.3.4 Effect of Root Thickness-Chord Ratio 

The effect of root thickness-chord ratio on wing stalling 
characteristics can be depicted from Figures 26, 27 and 28. These 
figures are presented for an aspect ratio of 6 wing having the 
tip thickness fixed at a constant value of 12% and the root 
thickness varying between 12% and 21%. Three families of airfoil 
sections are considered; 230 series, 44 series and 642 series. 

These figures show that for the wing employing 230 series 
airfoils increasing root thickness-chord ratio results in an 
increase in stall margins for outboard wing sections, a general 
inboard movement of the stall boundaries and a decrease in wing 
maximum lift coefficient. 

Similar effects of the root thickness-ratio can be noted for 
the wing employing 44-series airfoils, however, in this case 
the variation in the stalling characteristics is less pronounced. 

For the wing employing 642 series airfoil sections, an 
increase in root thickness-chord ratio causes a reduction in the 
stall margins on the outboard portion of the wing and tends to 
move the stall boundaries outboard in contrast to the trends 
indicated by the results for the other series. This effect is 
primarily due to the variation of c / max with thickness-chord ratio 
for this particular airfoil series. As can be seen from Figure 
20 the amount by which c/max decreases with increases in section 
thickness-chord ratios above 12% is greatest for the 230 series 
airfoils, less for the 44 series and least for the 642 series 
sections. In fact, for Reynolds numbers greater than 6x10° the 
642 series shows a small rise in c/max when the section thickness 
chord ratio increases from .12 to .15. 

Furthermore, for a given spanwise distribution of section 
Reynolds number and thickness-chord ratio the curve of c /max 
versus spanwise distance will be flattest for the wing employing 
642 series sections. This has the effect of promoting stall 
further outboard than would be the case for the wings composed 
of either 44 or 230 series sections. 

The fact that the dependence of c /max on (t/c) is least 
for the 642 series airfoils is again reflected in the results 
for wing maximum lift coefficient as shown in Figure 28. For 
this case increasing root thickness-chord ratio results in the 
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Figure 28. Concluded 
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smallest reduction in the value of maximum wing lift coefficient. 

5.3.5 Effect of Tip Thickness-Chord Ratio 

Figures 29 and 30 show the effects of tip thickness -chord 
ratio on wing stalling characteristics. The results are presented 
for a wing having a root thickness -chord ratio fixed at a constant 
value of 0.18 and tip thickness-chord ratios of 0.12 and 0.15. 
Limited computations were performed for tip thickness-chord ratio 
of 0.18 primarily for the purpose of establishing extrapolation 
trends for the variation of wing maximum lift coefficient. 

Figures 29 and 30 indicate that the effect of increasing the 
tip thickness-chord ratio is to reduce the stall margins and the 
value of maximum wing lift coefficient. This is due to the fact 
that an increase in tip thickness-chord ratio results in reduced 
values of section maximum lift coefficient in the vicinity of 
the wing tip, thus yielding lower stall margins. Also, the 
reduction in the values of maximum wing lift coefficient is 
smaller than that due to increasing root thickness-chord ratio. 
This is because the major effect of changes in tip thickness is 
confined to the tip regions where the wing loading is least. 

Although the effect of tip thickness-chord ratio on stall 
boundaries is not presented, it is expected to be similar to 
that of the root thickness -chord ratio discussed in the previous 
sub-section. 

5.3.6 Effect of Flight Reynolds Number 

Since the maximum lift characteristics of most airfoil 
sections are influenced by Reynolds number (see Figure 20) this 
parameter potentially represents a n effective means of controlling 
overall wing stall behavior. 

The effect of flight Reynolds number on stall margin distrib- 
utions, stall boundaries and wing maximum lift coefficient is 
shown in Figures 31, 32 and 33, respectively. These figures 
indicate that for wings employing 230 and 44 series sections an 
increase in flight Reynolds number results in increased stall 
margins over the outboard wing station and a shift of stall bound- 
aries towards the wing root. However, exactly opposite trends 
are indicated for wings employing 642 airfoil sections. An 
explanation of this behavior can again be obtained using the 
results of Figure 20. 

It can be noted from this figure that for the 230 and 44 ser- 
ies airfoils the rate of increase of maximum lift coefficient with 

Reynolds number is larger for the tip sections (t/c = 0.12) than 
that for the root sections (t/c = 0.18). This produces larger 
stall margins at the outboard (thinner) wing sections with increase 
in flight Reynolds number. 
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Figure 30. - Effect of Tip Thickness- Chord Ratio on C^gx 
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Figure 30. Concluded 
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However, for the 642 series sections Figure 20 indicates the 
opposite effect, namely that the rate of increase of maximum 
section lift coefficient with Reynolds number is much larger for 
the inboard sections (t/c = 0.18) than that for the outboard 
sections (t/c = 0.12). This produces a spanwise variation of 
maximum lift coefficient which has an increasingly downward slope 
toward the wing tips and, since the value of local lift coeffic- 
ient at any wing section increases with Reynolds number, the 
stall margins at the outboard sections will be reduced. 

The effect of Reynolds number on wing maximum lift coeffic- 
ient shows the expected increase with increasing Reynolds number. 
For the 230 series wing the maximum lift coefficient increases 
almost linearly with Reynolds number between 3 and 6 million, but 
with a further increase to 9 million the rate of increase is 
reduced. This trend is predictable from the section character- 
istics where it can be seen that the greatest changes in section 
maximum lift coefficient occurs below Reynolds number of 6x10^. 
For the wing series 44 and 642 the variation of maximum wing lift 
coefficient with Reynolds number in the range of 6 to 9 million 
is more linear, again reflecting the trend of the airfoil 
section characteristics 

5.3.7 Effect of Wing Camber 

The effect of a linear root-to-tip increase in wing camber 
on the stall margins is shown in Figure 34 for three values of 
taper ratio, three Reynolds numbers and zero wing washout. The 
results are compared with those for a constant camber. The 
camber variation chosen was 64218 root section and 64412 tip sec- 
tion. 


As can be seen from this figure the particular camber varia- 
tion chosen is not very effective in changing the stall margins 
on the outboard sections of the wing. This ineffectiveness is 
attributed, in part, to the particular choice of a linear varia- 
tion of camber. It is expected that a larger increase in stall 
margins would be obtained for a different combination of camber 
and wing airfoil sections. 

5.3.8 Effect of Fuselage 


Figure 35 presents the computer results for a wing alone and 
a high wing mounted on a fuselage of elliptical cross-section. 

The wing has an aspect ratio of 6 (based on the gross wing area) 
with zero washout, and root and tip thickness ratios of 0.18 and 
0.12, respectively. 

The computations were performed for taper ratios of 0.5, 

0.75 and 1.0 and wing-fuselage settings of 0°, 2° and 4°. The 
spanwise distributions of lift coefficient at stall were found 
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— Wing alone 

— Wing with Fuselage 
Washout = 0 Deg. 


Root Section 4418 
Tip Section 4412 
Aspect Ratio = 6 
Re = 6 x 1 Q 6 


Figure 35. - Effect of Fuselage 




to be nearly identical for all the wing-fuselage settings 
considered, thus resulting in the value of maximum wing lift 
coefficient ( C L max) being unaffected by the wing- fuselage 
incidence. The effect of positive wing- fuselage setting is 
merely to reduce the body angle of attack at which the stall 
first occurs. This reduction is approximately equal to the 
wing - f u se 1 age inc idence . 

Figure 35 indicates that for the high wing configuration 
shown, the effect of fuselage is to slightly reduce the wing 
loading specially in the region close to the fuselage. A similar 
reduction is indicated by the experimental results of Reference 
50 where the lift distribution on a high wing circular fuselage 
combination is presented. 

5.3.9 Effect of Partial Span Flap Deflection 

Figure 36 shows the effects of flap deflection and flap 
span on the stalling characteristics of a rectangular wing with 
64- series sections operating at a flight Reynolds number of 6 
million. The results are presented for the wing having a 20% 
chord split flap, deflected 60° and extending over 45%, 60% and 
70% of the wing span. 

It can be noted from this figure that, the deflection of a 
part- span flap lowers the stall margins on the outboard portion 
of the wing and causes the wing to stall further outboard. 
Furthermore, for the range considered, increasing flap span 
moves the stall boundaries inboard and increases wing maximum 
lift coefficient. It should be noted that the stall point on a 
flapped wing does not always occur at the flap end as would be 
predicted by simpler analytical methods than the one used in 

this program. 

A discussion of the above results and their influence on 
the design of an airplane for good stall characteristics is 
given in Section 7. 
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SECTION 6 


SCALE MODEL WIND TUNNEL TESTING 


It has previously been shown that the stalling character- 
istics of an isolated wing, with or without deflected flaps, can 
be adequately predicted through application of existing theoret- 
ical methods. In the case of a complete airplane, however, the 
interference of the fuselage, engine nacelles, propeller slip- 
stream, etc., on the flow over the wing may be such as to dras- 
tically modify the wing stalling characteristics. The interfer- 
ence effects of the fuselage and nacelles are predictable with a 
fair degree of reliability as long as potential flow conditions 
prevail. Unfortunately, the body interference effects often 
precipitate flow separation and available theories are not capable 
of predicting such phenomena. At the present time no theory is 
available to adequately predict the effect of the propeller slip- 
stream on wing stalling characteristics. 

For these reasons, it is highly desirable to obtain exper- 
imental information concerning the complete airplane stalling 
characteristics before the airplane goes into production. In 
certain cases economic considerations may indicate the construc- 
tion of a prototype model of the airplane with subsequent flight 
testing to obtain experimental information upon which to predicate 
the final design. In other cases, it may be more feasible to 
conduct a scale model wind tunnel investigation early in the 
airplane design stage in order to obtain the desired information. 

During the period prior to World War II wind tunnel tests 
were ordinarily made of models without propellers and empirical 
methods were relied on to account for the effects of propelller 
operation on the observed characteristics. Such a procedure was 
shown to be inadequate when quantitative flight test data became 
available. In consequence it is now considered essential that 
wind tunnel model evaluation of airplane flying qualities, 
whether they be concerned with stalling or with stability and 
control, should involve the use of a powered model. 

The following discussion will therefore consider some of the 
factors involved in powered model wind tunnel testing. 

6.1 SCALE MODEL REQUIREMENTS 


The selection of the model scale will be dependent on the 
size of the wind tunnel to be utilized in the investigation. As 
a rough rule of thumb the scale should be chosen such that the 
model wing span does not exceed approximately 75% of the width, 
or diameter, of the wind tunnel test section. Larger values of 
model size result in excessively large values of wind tunnel 
boundary corrections (Reference 23). 
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It is essential that the model be truly geometrically similar 
to the full scale airplane. In this regard, the leading-edge 
portion of the wing is particularly sensitive to deviations from 
contour. Caution must therefore be exercised to ensure that the 
airfoil shape over the forward 10 or 15 percent of the chord are 
true to the theoretical ordinates. Small deviations from true 
contour are not particularly critical over the aft portions of 
the airfoil. Wood, metal, plastic or combinations thereof may 
be satisfactorily utilized as material for model construction. 
However, it should be noted that composite wood and metal surfaces 
are not satisfactory because the wood shrinks and swells with 
change in atmospheric humidity conditions and thus gives rise to 
undesirable surface discontinuities. If wood construction is to 
be utilized, specially selected mahogany from the mainland of 
tropical America is suggested as being the most satisfactory. It 
has been found that mahogany from the islands of tropical America 
is not suitable for model construction because it shrinks and 
swells more and has a greater tendency to warp than mahogany 
from the mainland of tropical America. Similarly "Phillipine" 
mahogany is not suitable for model construction. 

6.2 PROPELLER DRIVE SYSTEM 

Specially constructed, compact squirrel cage induction 
motors have been utilized extensively as the propeller drive in 
powered models. Such motors are usually water cooled in order 
to increase their power output rating. In recent years some 
preference has been given to the use of compact pneumatic motors. 
The selection of motor type will depend to a large extent on 
availability of appropriate electric power supply or compressed 
air supply at the wind tunnel facilities being utilized. 

Propeller rotational speed can be measured through the use 
of a high precision tachometer and the measurement of propeller 
torque may be accomplished using an appropriate strain gage 
balance system. If a squirrel-cage induction motor is utilized 
it is possible to obtain a straight line calibration of torque 
versus minimum current. The minimum current point is obtained 
by varying the voltage-frequency ratio of the electrical power 
supply until the minimum current point is arrived at. 

The value of dynamic pressure at which power-on tests of an 
airplane model can be conducted is largely dictated by the torque 
rating of the propeller drive motor. For this reason it is 
advisable to at least make an approximate estimate of the critical 
simulation requirements before selecting the drive motor. 

6.3 SIMULATION OF POWER CONDITIONS 


In order to adequately simulate the effects of power in the 
wind-tunnel testing of models, it is essential that the axial 
and rotational velocity contributions of the propeller be in the 
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same ratio to the free stream velocity as that which prevails 
under the free flight conditions being simulated. This requires 
that the values of the propeller operating thrust and torque 
coefficients will change as the airplane flight speed or operating 
value of the lift coefficient is altered. It is usually most 
convenient to investigate the wind tunnel model through its 
operating range of lift coefficients at a fixed value of dynamic 
pressure, (constant velocity, fixed Reynolds number). To simulate 
the flight operating conditions during the constant velocity tests 
in the wind tunnel it is therefore necessary to provide a differ- 
ent operating condition of the propeller at each different value 
of lift coefficient or angle of attack investigated in the wind 
tunnel. 

Various techniques have been developed for satisfying the 
conditions of power similitude in wind tunnel testing. It appears 
that each different wind tunnel staff has its own preference as to 
the particular technique to employ. It is suggested, however, 
that the power matching technique described in Reference 21. 
is most appropriate for use in wing stalling investigations 
because that technique ensures a nearly exact condition of power 
similitude at each test condition and hence no interpolation of 
observed results is required. 

6 . 4 F LOW VISUALIZATION 

Numerous methods of visualizing the flow over the airplane 
model in the wind tunnel have been utilized. The most familiar 
are probably the smoke flow technique, some variation of the 
lamp black and kerosene method and the use of tufts. Of the 
various techniques that have been developed, tufting is the least 
complex and is usually the most satisfactory. In utilizing this 
technique, numerous tufts are attached over the upper wing and 
fuselage surfaces by cellulose tape or by other means. The tufts 
should be of flexible material such as wool or nylon yarn. The 
length of the tufts is not critical. Usually a tuft length of 
approximately 3 or 4 percent of the wing chord will be found 
appropriate. Tufts should not be located forward of 20% chord. 

The nature of the stall can be determined by noting the 
behavior of the tufts during the test conditions. Violent 
fluctuations and reversal of the flow direction as indicated by 
the tufts provides evidence of separation of the airflow from 
the surface under observation. The behavior of the tufts should 
of course be observed through a range of angle of attack from well 
below to well beyond the angle for maximum lift. 

The flow condition as indicated by the tufts may be recorded 
photographically, using either a still camera or a movie camera, 
or it may be recorded manually on the basis of visual observation. 
Each method has its own particular advantages. Attention is drawn, 
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however, to the fact that still photographs of the tuft flow 
pattern can be misleading. This stems from the fact that in 
the case of some configurations the nature of the airflow as 
the stall is approached may be very unstable and erratic. At 
one point in time, the tuft pattern may indicate the flow to be 
attached to the surface, an instant later the tuft pattern may 
indicate large areas of separated flow. If a still picture were 
obtained at the instant of attached flow it could lead to an 
erroneous conclusion. 

6.5 REYNOLDS NUMBER CONSIDERATIONS 


Stalling behavior of an aircraft can not be reliably predic- 
ted using small scale models in the wind tunnel. This is primarily 
due to the fact that it is extremely difficult if not impossible 
to duplicate in the wind tunnel the values of flight Reynolds 
numbers unless recourse is made to a compressed air tunnel or to 
a tunnel utilizing a high density gas such as Freon as a test 
medium. 

It has been well established that the maximum lift charac- 
teristics, including the wing stalling characteristics can be 
critically dependent on the value of the test Reynolds number. 

In consequence, judgment must be exercised in interpreting the 
wind tunnel stall test results in terms of the airplane flight 
Reynolds numbers condition. 

As an aid in interpreting the wind tunnel stall observations, 
it is suggested that the theoretical analysis described earlier 
in this report be applied to predict the stall at a Reynolds 
number corresponding to the wind tunnel test condition and at a 
Reynolds number corresponding to the airplane flight condition. 

By taking account of the differences between theory and experi- 
ment at the test Reynolds number, an improved estimate of the 
stalling characteristics at the flight value of Reynolds number 
may be obtained. This approach should at least give an indication 
as to whether the free flight stall condition will be more or less 
severe than the stall condition observed in the wind tunnel tests. 

6.6 MACH NUMBER CONSIDERATIONS 


The effect of Mach number on the maximum lift characteris- 
tics of airfoils has not been isolated’ and studied as thoroughly 
as the effect of Reynolds number. 

The results presented in Reference 51 indicate that at a 
given value of Reynolds number an increase in Mach number causes 
a moderate decrease in maximum lift coefficient under conditions 
such, that the local velocities on the surface of the wing are 
somewhat below sonic speed. It has also been clearly demonstrated 
by the results of References 9, 51 and 52, that when the free 
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stream Mach number is increased to the point that sonic speed is 
reached, locally on the wing, a large reduction in maximum lift 
coefficient takes place. The same results show that critical 
local velocities can occur at values of the free stream Mach 
number as low as 0.20. There have been instances in the past 
where wind tunnel investigators have attempted to improve the 
test value of the Reynolds number by conducting the tests at 
high values of the wind tunnel airspeed. Such a procedure can 
lead to quite misleading test results, particularly if the free 
stream Mach number is sufficiently high to permit the attainment 
of critical local velocities over the wing surface. 
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SECTION 7 


DESIGN PROCEDURES 

The results presented in Section 5 are intended to serve 
as a guide in the preliminary design phase of an unswept wing 
aircraft to determine the effects of wing geometric and aero- 
dynamic parameters on aircraft stalling behavior. While this 
data does not cover all the possible combinations of taper, 
twist, etc. which may be encountered, it should provide a basis 
for the assessment of the relative effectiveness of different 
wing designs in promoting acceptable airplane stalling character- 
istics. 

In the early design stage of an airplane the values of wing 
aspect ratio, taper ratio, and root and tip thickness ratios 
are usually chosen from considerations of performance, structural 
strength, etc. rather than stalling characteristics. In regard to 
wing performance, the computer program described in this report 
can be of value in providing data on wing lift, drag, and 
pitching moment characteristics through the complete angle of 
attack range, as well information on the span load distributions. 

The stalling characteristics of the basic wing design can 
be assessed from the design charts presented in Section 5, and 
if poor stalling behavior is indicated, the effectiveness of 
various methods for improving the stall can then be investigated. 
When a wing design emerges which promises to fulfill the perform- 
ance and stall requirements a final quantitative evaluation of 
its stalling characteristics can be made using the computer 
program which constitutes a part of this report. 

7.1 APPLICATION OF THE RESULTS OF THE PARAMETRIC STUDY 


The use of the design charts presented in Section 5 is best 
illustrated by a sample calculation described below. 

Consider a light, single engined airplane having the follow- 
ing characteristics: 


Wing aspect ratio 

Wing taper ratio 

Mean aerodynamic 
chord (m.a.c.) 

Wing loading 

Root airfoil section 

Tip airfoil section 


6 

0.5 

5.4 ft. 

4 

17.2 lb/ft ‘ 

23016.5 

23012 
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Cruise speed 
Cruise altitude 


154 m.p.h. 

10,000 ft., (standard day) 


A likely stalling speed for such an airplane, flans up, is 
about 70 m.p.h. or a Reynolds number, based on the m.a.c. of 
3.55 x 10°. A good estimate of the stalling speed for the 
untwisted wing is obtained as follows: 

(a) A number of speeds in the neighborhood of 70 m.p.h. is 
selected and the corresponding Reynolds numbers are calculated 
thus 


V ft/sec. = 

V m.p.h. = 

Re 


90, 95, 

61.35, 64.76, 

3.11, 3.28, 


100 , 110 , 

68.17, 74.98, 

3.46, 3.80, 


(b) Using Figure 33 for 0.5 taper ratio and a root thick- 
ness chord ratio of 0.18, the following values of c Lmax are 
obtained corresponding to the Reynolds numbers calculated in 
step (a), 


C L 


max t/c 



1.42, 1.435, 1.45, 1.47, 


(c) From Figure 28 the percentage change in c Lmax due to 
changing root thickness-chord ratio from .18 to .165 is estimated 
to be 2%. Strictly, the data applies only to Re = 6 x 10°, but 
the section characteristics, Figure 20, suggest that approximately 
the same changes can be expected at lower values of Reynolds 
number. 


( d) The increments in CL max obtained from step (c) are 
added to the c i-mox values from step (b) yielding: 


1.45, 


1.46, 


1.48, 


1.5 


max t/c = .165 

The results thus obtained are plotted versus stalling soeed 
in Figure 37. 

(e) A new series of values of c Lmax based on the relation- 
ship c Lmox = 2 x w/s are calculated thus: 
p v2 


Cl 


max = 


2 w/ s 

P v 


T = 1.78, 


1.60, 


1.445, 1 . 19 


These results are also plotted in Figure 37. 


(f) The intersection of these curves yields the values of 
Cl -max stall speed and Reynolds number as 1.475, 67.5 m.p.h. 
and 3.42 x 10°, respectively, as shown in Figure 37. 
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Figure 37. Variation of Wing Maximum Lift 
Coefficient with Stall Speed. 


L38 




At this value of Reynolds number Figures 31 and 32 show that 
the stall margin at the 70% station is only 0.02 and that the 
stall begins inboard at about 35% semispan and extends to about 
70% semispan. Obviously this is not acceptable and some means 
must be employed to move the stall inboard and increase the 
margin at the 70% semispan station. 

An acceptable stall pattern would be the one in which the 
outer edge of the stalled area began, say , inboard of the 40% 
station and the stall margin at 70% of the semispan was at least 
0.1. Figure 31 shows that, for the given taper ratio of 0.5, a 
stall margin of 0.1 at 70% semispan can be obtained with 7% 
degrees of washout. For this case the corresponding stall 
boundaries lie between 15% and 37% semispan. Since the use 
of washout influences the value of Cimax and hence stalling speed, 
steps (b) through (f) are repeated using the data presented for 
7% washout. If these computations are performed the new values 
of c Lmax , stalling speed and Reynolds number are 1.4, 69 m.p.h. 
and 3.5 million respectively. Using the new value of Reynolds 
number and washout of 7h° , Figures 31 and 32 indicate that 
satisfactory stalling characteristics of the selected wing are 
attained with the safe stall margin of 0.1 at 70% semispan and 
inboard stall boundaries extending between 15% and 36% semispan. 

If this amount of washout is used, a penalty may result in 
induced drag at the cruise speed. The magnitude of the increase 
in induced drag coefficient at the cruise lift coefficient of 
.38 and cruise Reynolds number of approximately 6 million is 
obtained from Figure 25 as A Coi= .0012. The significance of this 
drag increase depends on the drag coefficient of the complete 
airplane. For an airplane having a drag coefficient equal to 
0.02, an increase in induced drag coefficient of 0.0012 represents 
a 4 m.p.h. reduction in cruising speed at the cruise power setting. 

An alternative to using a large amount of washout as 
high as 7h degrees is to incorporate a tip section of 
higher camber than the root section with linear fairing in 
between, e.g. change the tip section from 23012 to 43012. On 
the basis of the results presented in Figure 34, for a wing 
employing linear camber increase from root to tip, it appears 
that increasing camber alone will not result in any significant 
improvement in the stalling characteristics. Combinations of 
linear camber increase with moderate amounts of washout may 
result in an effective compromise to yield satisfactory perform- 
ance and acceptable stalling characteristics of the selected wing 
configurations. The parametric investigations of such effects 
can be easily performed utilizing the computer program presented 
in this report. 
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If an effective combination of camber increase and washout 
cannot be found the only remaining wing parameters which might 
influence stalling characteristics are aspect ratio, thickness 
distribution and taper ratio. 

On the basis of the results shown in Figure 23, changing 
aspect ratio is ineffective. Increasing root thickness-chord 
ratio from .165 to say .18 is also of little value according to 
the results presented in Figure 27. 

While an increase in wing taper ratio would represent a 
major change if the wing design were sufficiently far advanced, 
it might be less expensive, in the long run to make such a 
change than to try to solve bad wing stall problems by other 
means during the flight testing phase of development. 

An illustration of the strong influence which is exerted 
by taper ratio on wing stalling characteristics is provided by 
repeating the above calculations for a taper ratio of 0.75 with 
2%° of washout. The increased taper ratio is achieved by reduc- 
ing the root chord and increasing the tip chord by the same 
amount so as to maintain the same wing area as in the original 
design. 

By interpolation between these results and results obtained 
for taper ratio of 0.5 it is found that an increase in taper 
ratio from 0.5 to 0.65 will produce a stall margin at 7/10th 
semispan equal to 0.1 with the stalled area extending between 
10% and 40% semispan. The calculations also show that the change 
in taper ratio and the incorporation of 2^° of washout does not 
alter the stalling speed to any significant degree. Furthermore, 
the induced drag penalty due to washout for this taper ratio 
is negligible. 

If all of the above measures fail to indicate acceptable 
stalling characteristics then it must be left until the flight 
test phase to try to improve matters by the use of the various 
"fixes" discussed in Section 2. Even at this stage the computer 
program should prove valuable in assessing the relative merits 
of the post design modifications. For example, if the install- 
ation of sharp wedges over a portion of the leading edge is 
being considered^ precise evaluation of the effectiveness of 
such a device can be made by using the computer program if data 
is available pertaining to the aerodynamic characteristics of 
sharp nosed sections. Some data on the effect of sharp leading 
edges on section maximum lift coefficient can be found in 
References 53 and 54. 

7.2 GENERAL CONSIDERATIONS 

The design procedures discussed in this section together 
with the stall charts presented in Section 5 are considered 
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adequate for most practical cases in evaluating stalling character- 
istics of unswept wing aircraft. These procedures and charts 
should be specially valuable in the preliminary design phase in 
which numerous trade-off studies are required to obtain the best 
compromise in aircraft configurations. Once a given design has 
been frozen ensuring satisfactory aircraft performance, stalling 
characteristics, handling qualities, etc. it is recommended that 
the computer program presented as part of this report be utilized 
to more accurately predict the stalling behavior of the final 
aircraft configuration. As a by-product, the computer program 
will also yield valuable performance information such as 
distributions and integrated values of wing lift, drag and pitching 
moment coefficients, for cruise or any other aircraft operating 
condition. 
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SECTION 8 


CONCLUSIONS AND RECOMMENDATIONS 

Using the results presented in this report, the following 
conclusions and recommendations are made: 

1. Based on good correlations obtained between the theo- 
retical results and the available test data, it is concluded that 
the lifting line theory can be confidently used to predict stall 
characteristics of wings having aspect ratios of 6 and larger. 
This theory is expected to yield satisfactory predictions of 
overall load characteristics for wings of aspect ratios as low 
as 4.0. 


2. From the results of the parametric study it can be 
concluded that taper ratio is one of the most effective design 
parameters influencing aircraft stall characteristics . Increase 
in taper ratio results in an increase of the stall margins on 
the outboard sections of the wing and in an inboard shift of 
stall boundaries. This, however, is accompanied by a reduction 
of maximum wing lift coefficient. 

3. Washout may be used to improve stalling characteristics 
of moderate to high taper ratio wings. This improvement, however, 
may entail a performance penalty associated with a reduction of 
wing maximum lift coefficient and an increase in wing induced 

and profile drag. 

4. An increase in root thickness -chord ratio and flight 
Reynolds number yields favorable effects on wing stalling 
characteristics associated with the 230 and 44 series airfoils. 

For the 642 series sections such increases result in unfavorable 
effects. Increasing tip thickness -chord ratio has an unfavorable 
effect for all three airfoil series. For all cases the values of 
wing maximum lift coefficient are reduced. 

5. For the wing configurations investigated in this 
report, deflection of part-span flaps lowers the stall margins 
on the outboard portion of the wing and shifts the stall bound- 
aries outboard. Increasing flap span moves the stall boundaries 
inboard and increases wing maximum lift coefficient. 

6. The effects of aspect ratio, linear camber, and fuselage 
on wing stalling characteristics appear to be small and may be 
neglected for most present-day light aircraft. 

7. The analysis and the design charts presented in this 
report apply strictly to unpowered flight, out of ground effect, 
as would be the case in aircraft approach to landing. 
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8. Based on the work accomplished in this program, it is 
recommended that the theoretical analysis be extended to include 
propeller slipstream and ground effects. Furthermore, additional 
flight test and wind tunnel test data should be obtained for the 
purpose of verifying the theory. 
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APPENDIX A 


INTERNAL LISTING OF THE COMPUTER PROGRAM 


0VCftLAY(lJLIMDA,0,-O) 

PROGRAM S FALL (. I NPUT = 20 L , OUTPU T = 100 1 , TAPE 8 = 1 NPUT » TAPE6=0UTPUT , TAPE 1 

1 = 201 , I' APE 2 = 201 , TAPE 3 = 201 , TAPE 4 = 201, -TAPE5 = 201 , TAPE 1 0 = 20 1 , T A PE 1 5 = 20 1 

2 , T A PL 20 = 20 1 , T m P 244=20 1 ,'T APE 7 = 1 00 1 « TAPE 9 = 2 U 1 ) 

PROGRAM- STALL COMPUTES THE STALL CHARACTERISTICS OF STRAIGHT WING 
AIRCRAFT 

I HIS PROGRAM WAS PREPARED FUR NASA LANGLFY RESEARCH CENTER UNDER 
CONTRACT NA S 1 — 8 380 ttY UYNASCIENCES CORPORATION, SCIENTIFIC SYSTEMS 
DIVISION,' BLUE KELL, PENNSYLVANIA 

THE DYN A SCIENCES PERSONNEL WHO CONDUCTED THE INVESTIGATION WERE 
MESSRS. M. A. MCVEIGH, II . KISIELOWSKI, AND J. G. MCHUGH. AREA CODE 
2 15-65 3-0250. 

ANY ERkOnS Ok PROBLEMS ENCOUNTERED IN USING THE PROGRAM SHOULD BE 
DIRECTED TO MR. CHARLES H. FOX, JR. AT NASA LANGLEY. AREA CODE 703 
-827-3711. 

A CARD DECK AND DOCUMENTATION FUR THE PROGRAM ARE AVAILABLE FROM 

cosmic, university of gegkCia, Athens, gforgia, 30601. 

THIS Put GRAM- IS WRITTEN IN CDC FORTRAN IV, VERSION 2.3, TO RUN ON 
CDC 6600 SERIES COMPUTED WITH THE SCOPE 3..0 OPERATING SYSTEM AND 
LIoRARY TAPE.. THIS PROGRAM WAS WRITTEN TO PERMIT CONVENIENT 
MODIFICATION TOR USE UN OTHER COMPUTERS. 

OIMrNSIl-N ARRAY! 5,25, 12 ) ,C(.19 ) ,EPS( 19) ,'TKANS( 19 ) ,REY< 19) ,ETA( 19 ) , 
1H0PPI L9) ,CLMAX M9 ) , ZHERL ( 2 ,'6 ) , WHERE! 2,6) ,Y( 19 ) , TAU( 19) ,BETA ( 19, 19) 

2, TRIX! L9 , 19 ) ,M.AZZ (.6 ) ,MZCOL (.6) ,-MAXX (.6 ) , MXCOL < 6 ) , XHERE ( 2 , 6 ) , M'wCOL ( 6 ) 

3, MAw*(-6 ) ,CM( 19) , CBG! 19) ,CVAL ( 19) , ALPCI 19) ,CBC< 19 ) , ALPIIU! 19 ) , ALPH! 
419 ) , ALPHZ ( 19) ,'ALPHE ( 19 ) ', DELTA (.19 ) ,CL ( 19 ) , ALPHV! 20) , YHERE (2 ,6 ) ,MYCO 
5L (.6 ) ,MAYY ( 6 ) ,CL ADD (.19) ,CLDFL ( 19 ) , CLA02 (19) , CLAD! ( 19 ) , F ( 19 ) , ALPC ( 19 
6) ,■ TON Y ( 19) 

COMMON I NNOW , ISwI T (.3) , ALPHA ,kLYN, CL L ,-RE YON, XMAX , ALMAX ( 19 ) , CLMAX ,C , 
1EPS,TkA,\S,REY,'ETA,H0PP,ZHERF,WTiERE, aRRAY, Y, BETA, TFAC, TRIX, TAU,MAXX 
2,M'AZZ,MXC0L,MZC0L,aSPEC, rAPbRjBF, REYNOLD I SC R, PIER, CRB, Q,TSTAX»tDGL 

3 ,'S I G , aI PHR ,'nF I. A P , NL VL ,'NP , IY, IZ, IR, I P ,-I S I S , I STAR , A , B , H , TAUT , TAUR , 

4 TWIST , R, R W X, YF'E RE, M YCCL.MAYY, FLAP, TON Y ,.T W I S A , X , Z , CM, ACC, XHERE, MWCO 

5 L , MAUL , C A MB (19) , CAM f.R , CAM.bT , DUMY I , D0MY2 , MflMt ( 25 ) , AHERC (2,6) , MA A A ( 6 
(,) » M; A C L L ( 6 ) ,'BHERE (2 , 6 ) ,'MAhli (.6) ,MBCUL(6) ,'CHERE(2,6) ,MACC(6) ,MCCC)L(6) 
7 .DHEkE ( 2 ,o ) ,M.ADLM 6 ) , MDCUL (.6 ) , ALPHF, 'ALPHV , CL ADD , C LDE L , C L AD2 ,CLAD1 , 

F , 1 A I ,FF ,LDCER, IPRAB, IKT, IOT , l AGIN, 1 TR, I S ,L;UM< 100 ) 

BL I NUA = 6LBL I NO A 
RECALL=6hRECALL 

1.0 CALL OVERLAYItJLINDA, 1 ,'U , RECALL) 

IF(.MLAP): 20,20,30 
20 CALL OVERLAY ( BL INL:A, 2 , 0, RECALL) 

GU TU 10 

30 CALL OVFRI.AYIBLINDA, 3,0, RECALL) 

GO T (j 10 
END 

SUBROUTINE M INV ( A ,-N , D ,'.L , M ) 

MATRIX INVERSION SUBROUTINE 
DIMENSION A(1),L(1>,M(L) 

SEARCH FUR LARGEST ELEMENT 

0=1 . 0 
NK = -N 

DO 160 K - 1 , N 

NK=NK+M 

L(..K)=K 

M(-K)=K 

KK=NK+K 

GIGA=A(KK) 
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DO 20 J = K,N 
I Z=N* ( J— 1 ). 

DO 20 I=K,N 
IJ=IZ+I 

IF ( ABS ( B IGA )-ABS ( A (• I J ) ).') 10,20,20 
1.0 B I GA=A UJ) 

L IK ) = I 
M (.K ) = J 
20 CONTINUE 

INTERCHANGE ROWS 

j=*uk) 

IF(J-K) 50, '50, 30 
30 K I = K— N 

DO AO 1 = 1, N 
KI-KI+N 
H0L L)=-A ( K I- ) 

J I * K I — K + J 
A|,KI)-A(JI) 

A0 A (■ J I )=H0L0 

INTERCHANGE COLUMNS 

50 I ■= M ( K ) 

IFM-N) 80, BO, 60 
60 JP.= N*H-1 ) 

DO 70 J = 1 , N 
JK=NK+ J 
JI=JP+J 
H0LU=-A( JK) 

A(.JK)=A( JI ) 

70 A(.JI)=H0LU 

DIVIDE COLUMN BY MINUS PIVOT (VALUE OF 
PIVOT ELENiENT IS CONTAINED IN BIGA) 

80 IF(ABS(BIGA)-l.t-20) 90,90,100 
90 D=0 , 0 
RETURN 

100 00 120 I - 1 ,'N 

IF(I-K) 1L0, 120, 110 
1L0 I K = N K + I 

A (.IK) =A( IK) /(-BIGA). 

120 CIJN 1 INUL 

RtDUCc MATRIX 

UO 150 I — L , N 
1 K=NK+ 1 
H(ULD = A (IK) 

I J= I— iN 

DO 150 J=L,N 
I J=I J+N 

IF(I-K) 130,150,130 
130 IF(.J-K) 1 AO ,150, 1 A 0 
1A0 K J= I J- I+K 

A (. 1 J ) = H U L D * A ( K J ) + A (. 1 J ) 

150 CONTINUE 

DIVIDE ROW BY PIVOT 

K J = K - N 

DO 170 J = 1,'M 
KJ=KJ+N 

If(.J-K) 160,170,160 
"160 ATK J )'= A ( K J 7 / KI G A’ " 

170 CONIINUF 
C 
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PRODUCT OF PIVOTS 


U=D*b 1 GA 


REPLACE PIVOT BY RECIPROCAL 

A (.K K ) = 1 * 0 / B I G A 
180 CONTINUE 

FINAL ROW AND COLUMN INTERCHANGE 

K = N 

190 K-K-l 

I F ( . K ). 260,260,200 
200 I =L ( K ) 

IF( I-K). 230,230,210 
2 L0 JQ=N* { K— 1 ) 

JR = N*M-1) 

DO 220 J = 1,N 
JK= JQ+ J 
H0LD=A(.JK) 

J I = JR+ J 
AUK)=-A(JI ) 

220 A (. J I ) = W0LD 
230 J = M ( K ) 

I F ( J — K ) 190, 190,240 

240 K I = k-N 

DU 250 1 = 1, -N 
K I - K I + N 
HDLD=A(KI ) 

J I =k i-k+. J 
A {.K I )=-AUI ) 

250 A (. J 1 ) = HO L D 
GO TU 190 
260 RETURN 
END 

SUBKUU1 INF DAGET( ARRAY, I FILE, I I ) 

SUBROUTINE TO GET TABLE FROM DISK AND PUT 
IT INTO CORE 

D I MUMS ION ARRAY ( 5 ,2*5, 12 ) , D R R A Y. ( 5,25,12). 

COMMON I MNOW, I Sh I I' (.3 ) , ALPHA, REYN, CLL,«Rfc YON ,«XMAX f ALM AX ( 19) , CLMAX ( 19 
1 ) , C ( 19 ) , E PS { 19 ) ,T RAMS ( 19 ) ,-Rfc Y (.1 9 ) , E TA ( 19 ) ,HOPP ( 1 9 ) , ZHERL ( 2 ,*6 ) , 
2wHE<E(2,6) ,-DRRAY, Y(.l9) ,BETA( 1 9 , 1 9 ) ,TF AC , TR I X ( 1 9 , 1 9 ) ,TAU( 19) ,MAXX(6 

3 ) ,*MAZ2(.6) , MXC0LI6 ) ,’FZ COL ( 6 ). , AS PEC , T APFR , BF ,.R E YND , C I SCR , P I ER , CRB , Q , 

4 rSTAX , EDGE ,«S I G , AL PHR , NF L AP , (ML VL , NP ,■ I Y. , IZ,IR,IP,ISIS,ISrAR,A,B,H, 

5 TAUT, IAUR, TWIST , R , B X , Y H t R E ( 2,6) ,-hYCUL (.6 ) ,MAYY(6),FLAP,T0NY(19), 
6TWISA,X,Z,CM( 19) , ACC ,‘XHE R E (■ 2 , 6 ) , NWC0L(6) ,MAWW(6) ,CAMB( 19) , CAMP R, 
7CAMBT, DUM Y1 , DUM Y2 ,-NAME (.? 5 ) A 1 1 1 K E ( 2,6 ).,MAAA{ 6 ) , MACOL (6 ) 

TABLE PRESFNTLY IN CORE IS INNOW 

IF TABLE IS ALREADY IN CORE THEN RETURN 
IF NOT THEN GET TABLE FROM DISK 

I F ( IF ILL- INNOw) 10,90,10 
10 11=1 

REWIND KILL 
RE AD ( IF ILE ) ARRAY 
REWIND IFILE 

RESET TABLE PRESENTLY IN CORE INDICATOR 

~i • 

I N N D W - I F ILE 
IF CNF LAP) 90,20,90 
20 IFMFILL-5) 50,30,90 
30 DO AO J FUX- I , ML VL 
K R 0 w = M A A A ( J F 0 X ) - 2 
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KC 0 L - M A C GL (. J F 0 X ) 

DO AO KB0B=2,KCGL 
40 ARRAY ( JFGX, KRt;W,KbOE )=0. 

GO IU 90 > 

50 IF MFl’Lfc— 1) 90,60,90 
60 DG EC J F G X - 1 » N L V L 
KROW-MiAXX ( J FOX 1-2 
K C 0 L = f' X C CJ.L (.JFGX) 

DG 70 K B 0 1 * - 2 t K C U L 
70 ARRAYt JFCX,KRCW,KBGB )=0. 

00 CP NT 1 N.UE 
90 RETURN 
END 

SUBROUTINE AAA (.CL I S F , NV ) 

SUBROUTINE Til OUTPUT A ONE DIMENSIONAL 
ARRAY 

DIMENSION CLISUL) 

I P = 6 

WRITE ( IP, 10) (.J, CL 1ST! J ) ,J=1,NV) 

WRI-TU IP, 20) 

RE TURN 

10 FORMAT ( 5 (. L X , 1 H (. , I 2 , 2 H ) ,F16.8>> 

2 0 FORMAT ( 2 UH ) ) 

END 

SUBROUTINE Z Z 7 (.VALUl ) 

SUBROUTINE TU OUTPUT A SINGLE VALUE 
OUTPUT VALUE IN bUTH F AND E FORMAT 

IP = 6 

WRI-Tt ( I P , 10 ) VALUE ,< VALUE 
RETURN 

10 FORMA T (. IX, FI 5. 8, El 6. 8/ IX ) 

END 

SUBROUTINE SSSUtWP) 

SUBROUTINE TO OUTPUT A SQUARE TWO 
DIMENSIONAL ARRAY 

DIMENSION A ( 3 9, L 9 ) 

I P = 6 

WRI TE ( IP, 10) (. {, J , K , A ( J,K ) ,K=1 , N P ) , J = 1 , N P ) 

W R I T E ( IP, 20) 

RETURN 

1.0 FORMAT (.5(.LX, lH(. t 12, 1H, I2,2H) ,E16.8) ) 

20 FORMAT ( 3 (■ LHO/ ) ) 

END 

FUNCTION FEKP <.R1 ,«R2,R3, Cl ,'C3 ) 

INTERPOLATION FUNCTION 

I F ( R 1 — R 3 ). 20,10,20 
10 TE R P =C 1 
RETURN 

20 TE RP=£ 1+ ( (-C 3-C L ) *. ( CR2-RL J/IR3-R1 ) ) ) 

RETURN 

END 

SUBROUTINE StTSW 


SUBROUTINE CHECKS SWITCH SETTINGS UN DATA 
CARD AND LISTS OPTIONS SELECTED 


COMMON INNOW, I SW I I (.3) 
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K=0 
I P=6 

w RME( IP, 90) 

I F ( ISWIT(.L) ) 20,20,-10 
10 K=K+ 1 

wki rt( ip, loo) 

20 IFUSU-H2) ) 40,40,30 
*0 K = K+ 1 

W R I* 1 t ( IP, 110) 

40 IFUSWT-l (.3) ) 60,60,50 
SO K=K+1 

WRITE! I P, L20) 

60 if(K) 80,70,60 
70 WkITCMp, 1 3 r * ) 

80 RETURN 

90 FORK AT (. 1H 1/ 1H0 / 1H0 ,’50X , 29HSP tC I AL FEATURES FOR THUS RUN/1H0) 

100 FORK AT (.SIX, 14HRE AD IN VALUES/1H0) 

1L0 FORMAT {.Six, 34HTRACE (OUTPUT OF KEY C A LCUL AT 1 ONS ) / 1H0 ) 

120 F0RMAT(.51X, 28HPRINT MATRICES SETA AND TRIX/iHO) 

130 FORMAT TblX, 31FN0 SWITCH SETTINGS REGULAR RUN) 

b N D 

SUBROUTINE DAT Srt < N , I CGND ) 

SUBROUTINE MAKES DECISION BASED ON THE 
SETTING OF THE SWITCH SLT=1,N0T SET=0 
I C 0 N D = 1 IF SWITCH IS SET 
I COND = 2 IF SWITCH IS NOT SET 

COMMON INNOW, ISwIT ( 3) 

IF C.N-3 ) 10,20,30 

10 IF ( IS* IT (.1) ) 40,40, 50 
20 IFI.ISaI 1 (.2) ) 40,40,5 0 
30 IF ( ISWI T (.3) ) 40,40,50 
40 ICGND=? 

R b T U K N 
5 0 l C 0 N L) = 1 
RE TURN 
END 

SUBROUTINE AF.RD A (.ARRAY , NLVL , WHERE , MXCOL , MAXX , I FILE, II ) 

SUBROUTINE TU READ AERODYNAMIC TABLES FROM 
CARDS AND PUT TABLES ON DISK 

DIMENSION NAME (.40 ) , MAXX ( 6 ) , MXCUL < 6 )!, WHERE ( 2,6),ARRAY(5,25,12) 
COMMON 1 IMN Ow 

M A X X ( 6 ) CONTAINS NUMBER OF ROWS IN A TABLE 
W H E R F ( 2 , 6 ) CONTAINS TABLE NUMBER AND TAU 
VALUE 

INNOw IS THE TABLE PRESENTLY I-N CORE 

I R = b 
1P = 6 

INITIALIZE LEVELS IN WHERE 

DO 10 J=1 ,NLVL 
1.0 WHFKE( L , J ) = J 

LOOP TO READ IN LEVELS OF A SET OF EITHER 
LIFT, DRAG, OR PITCHING MOMENT 

DO 40 L V L = 1 ,NI.VL 

READ AND PRINT NUMBER OF ROWS »• COLUMNS, TAU 
VALUES FOR A GIVEN LEVEL 

READ ( I R, 50) NC,NCOL, WHERE(.2,LVL ) 
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c 

c 

c 

0 

0 

c 

c 

c 

c 


c 

c 

c 


c 

c 

c 


c 


c 

c 

c 


c 

c 

c 

c 


c 

0 

c 

c 

c 


WRITE <- I-P ,60 ) NC, NCOL, ’WHERE (.2, LVL ) 

STORE NUMBER OF COLUMNS 

MXCUL(LVL)=NCOL 

STORE NUMBER UF ROWS 

M A X X { 1 V L ) =NC 

READ AND WRITE TITLE OF TABLE 

READ UR, 70) NAME 
WRITER I P ,80 ) NAKfct-LVL 

READ VALUES FOR TABLE 

DO 20 1=1, NC 

20 READ( I R , 9 0 ) (ARRAY (.LVL , I , J > , J = 1 , 1 2 ) 

DO 30 J J “ 1 , NC 

30 in R I f E UP, 100) (-ARRAY { LVL , JJ, I I ) , I 1 = 1 ,<NCuL ) 

AO CONTINUE 
I 1 = 1 

W Rift COMPLETE TABLE UN DISK 

WRITE ( I E ILE ) ARRAY 
RETURN 

100 FO Rt AT (■ L2FI0. B) 

50 FORMAT (212, El 0.0) 

60 F 0 R M A T (. 1 X , 2 I 2 , F 5 . 2 T 
70 FORMAT LA 0 A? ) 

80 FORMAT (.1 1X,40A2,2X, 12) 

93 FGRMAT(.E7.3,9F8.3> 

END 

SUBRUU1 INL L> K I D G ( C C Z 1 , N S , N P , L D M , I L XY ,’I E , E'Y,*AU , AMB , XX, I S , I P ,YY ) 

SUBROUTINE TO INTERPOLATE BETWEEN TABLES 

DIMENSION XX { I ) , EYH ) , AU( I ) ,CCZ 1 ( 1 ) ,CCZ2 ( 19 ) ,CCLR1 ( 19),CALR1(I9), 
ICC LX 2 ( 19) , CALK2 ( 19 ) , A IMP (19) , R£Y( 19 ) , TAU( 19 ) , ARRAY ( 5,25, 12 ) ,C ( 19) , 
2LPS ( 19 ) , TRANS (.19) , EI AU9 ) ,HOPP( 19) ,.CLMAX (19 ) , ZHERE ( 2,6) , WHERE (2 ,6) 
3, YU9) , BE TA (19,19 ), I RIX (19, 19 ) ,MAZZ ( 6 ) ,MZCOL (6) ,MAXX (6) ,MXCOL ( 6 ) , 
4YHERE ( 2 , 6 ) , MYCOL (6 ) , MAYY ( 6 ) ,CM{ 19) , XMERE ( 2 ,6 ) ,MAWW( 6) ,MWCOL ( 6 ) 
COMMON INNOW, I S w I T (.3 ) , AL PHA , K l YN , CL L ,*RE YCN , XMAX , ALM AX ( 19) ,CLMAX,C, 
1EPS, TRANS, RirY.tT A, HOP P, ZHERE, WHERE,' ARRAY, Y, BETA, TFAC, TR IX, TAU,MAXX 

2 , M A Z Z , M X C 0 L , M7C0L ,-A S PEC , f APER , B F ,RE YND ,D I SC R, PIER, CRB, 0, TSTAX,EUGE 

3 ,‘S i-G, ALPFR,nFLAP,NLVL,NQ, IY, IZ, IR, IQ, IT, I STAR, A , B , H , T AUT , TAUR , 
4TWI-S1 ,R,BWX,YHL.<t,MYCOLf N AYY, FLAP, TONY U9) , T W I S A , X , Z , CM , ACC , XHERfc , 
SMWCUL, MAW W, -CAMP ( 19) , C A M B K , C A M b I , DUMY1, UUMY2, NAME ( 25 ), AHh RE ( 2 , 6 ) , 

6 M A A A ( 6 ) , M A C 0 L (. 6 ) , BULK L (-2 , 6 ) , M ABB ( 6 ) ,MBC0L(6) , CHfcRE < 2, 6 ) ,MACC(6) , 
7MCCUU6) ,UHLRF(.2,6) , M ADD ( 6 ). , MOGUL ( 6 ) 

INITIALIZATION 

IF R E Y 0 N = 9 9 9 THEN SPECIAL CASE 

DUi‘;Y2 = 0.0 
DUKY 1=0.0 

I f-< -K L Y UN- 9 9 9 . ): 20,10,20 
10 K E Y = 1 
CO TO 30 
20 KL Y = ? 

STORE C VALUE 

30 C L T = C L L 

NEW IS TF!E FILE NUMBER OF THE OTHER CUBE 
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PRIMARY CUBES ARE NUMBERED 1, 2,3,4 
SECONDARY CUBES ARE NUMBERED 5.10,15,20 

NFw -L0M*5 

STORE XMAX 

XMX=XMAX - -- 

DETERMINE IF SINGLE VALUE IS TO BE USED OR 
LIST OR VALUES IS TO BE USED IN LOOK UP 

IF( 15-IP) 50,40,50 
40 K G 0 = 1 

SET UP FOR CONSTANT VALUE OF X 

A T M P ( 1 ) = X X (. 1 ) 

GO HJ 70 

SET UP FOR VARIABLE VALUE OF X 

50 KGO = 2 

DO 60 J = N S , N P 
6 U ATMP ( J ) = XX <• J ) 

PUT PROPER TABLE IN CORE 
70 CALL OAGE M, ARRAY, NEw, IZXY ) 

SET UP ALPHA EITHER VARIABLE OK CONSTANT 

DO 1 9 U K = N S , UP 
GU TO (b 0,90) , KGO 
80 AL PHA = A T M;P (.IS) 

gu ru l.oo 
90 ALPHA=ATMP (.K ) 

GET TAU VALUE FROM ARRAY 

100 T A U X = T A U (■ K ) 

GO 1C ( 12 (if 11**) , KEY 

GET PRUPFR REYNOLDS— NUMBER EITHER REGULAR 
OR MAX NUMBER 

110 KEY ON = R t. Y ( K ) 

RE Y\=999 . 

GO ro L 30 
12 U RE YU = Kfc Y (.K ) 

130 C L L = C L r 
X M A X = X i v X 

PERFORM LOOK FOR LIFT, DRAG, PITCHING 
MOMENT OK ELAP CASE 

GO TO ( 14U, 15C, 160,1 70)!, LUM 
140 CALL A R C (. A R R A Y , T A U X , Mi A A A , M A C U L , I E , A H E K E » \) L V L ) 

GO re, 1.8 0 

15U CALL ARC (.ARRAY,- T AUX , M ABO , M BC CL , It, P HERE, NLVL ) 

GO 1(1 1>.Q 

160 CALL ARC LARRAY ,-T AU£ , MACC , MCCUL , I E , CHbRE , NLVL L 
GO lu Lt-O 

170 CALL AKC(AKRAY,TAUX,.VADD, MOGUL, 1 1 , UHf RE , NLVL ) 

STUKL VALUES FOUND IN LOOK UP 

180 CC Z 2 ( K ) = Y Y 

CCLR2 (K ) =DUN Y 1 
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C A L N 2 ( K ) = DU M Y 2 


E P li: A r FOR NUMBER Of VALUES TU BE FOUND 


190 CONTINUE 

GET Nt-XT TABLE IN CORE 
CALL OACETt ARRAY, LDK, IY) 

Sc T UP ALPHA VALUE 

DO 310 K=MS,NP 
GlJ TO (200,210), KGO 
?0C ALPDA = AIMP(. IS). 

GO TO- 220 
210 ALPHA = ATMP(.K) 

SET UP REYNOLDS NUMBFR 

220 GO TO (240,230), KEY 
230 RE'Y0N=KE Y ( K ) 

REYv=999. 

GO TO 2 ft 0 
240 REYN = NEY(.K) 

SET UP C V A L ,’ TAU, AND XMAX VALUE 

250 CLL=CLT 

TAUX = TAu (.!< ) 

X M A X = X M X 

LOOK UP FOR SECOND LIFT, BRAG, PITCHING 
MOMENT, OR FLAP CASE 

GO TO ' 260 ,'270, 280, 290 ) , LDM 
260 CALL ARC (.ARRAY , T AUX,KAXX,MXCOL, I E, WHERE, NLVL ) 

GO TU 300 

270 CALL ARC URR AY, TAuX.MAZZ ,MZCOL, I E, ZHERE , NLVL ) 

GO TO 300 

280 CALL ARC (.ARkaY , T AUX , MAYY , MYCOL , I E , Y.HERE , NLVL ) 

GO TO 300 

290 CALL ARC ( ARRA Y » T AUX , MAWN , MWCCL » I E * XHERE , NLVL ) 

STORE VALUES POUND IN LOOK UP 

300 CCZ 1 ( K ) = YY 

CCLR 1 ( K ) =UUi v Y 1 
CALRK.K )=0UMY2 
310 CUN I liiUL 

SPECIAL CASE FUR LIFT LOOK UP 
GO TO (32u, '400, 40 0,400) , LDM 

REPEAI FOR ALL VALUES OF TAU 

320 DC) 390 K=TMS,NP 
ALPHA = ATMP(.K ) 

CCLR = TEKP(CAM,BR,CAML1(K) ,CAMBT,CCLRL(K) ,’CCLR 2 ( K ) ) 

CALR= TE'RP (CAMHR.CAMB (K ) »CAMBT ,-CALR 1 IK) ,'C ALR2 ( K ) ) 

IF (.( ALPHA— CALRKK) )>( ALPHA-CALR2 ( K )! ) ). 340,330,330 
3 30 CCZ l ( K ) =TERP ( CANBK.CAMB ( K ) , CAMBT ,CCZ1 ( K ) ,CCZ2 ( K ) ) 

GO TO 3 ‘>0 

340 IF(CaLK2(.K)-CALk1(K) ) 350,350 ,'360 
350 I F (.ALPHA-CALK ) 370,370, 380 
360 I F ( ALPHA-CALK ) 380,380,370 

370 CLA=CU<P (.CALRl(.K) , CALR2 ( K ), ALPHA , CCLR 1 (.K ) ,CCZ1(K) ) 

C L C = 1 fcptP (.CAKhR » C AMb ( K ) ,CAMBT , CLA ,.CCLR2 (-K ) ) 
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CCZl(K) = Tf:RP(CALR?(.K ) , AL PHA , C ALR , C LC ,CC LR ) 

GO ru 100 

380 CLri = I'LRP (.CALR2 (.K ) , C ALk L ( K ) , AL PHA , CC LR2 (.K ) , CC Z 2 ( K ) ) 
C L C = T E R P (X A M [3 K » C A MB ( K ) ,-CAWfJ ( , CC L R 1 U ) , C L B ) 

CCZ 1 ( K ) = T E R P ( C A L R 1 (-K ), ALPHA , CAL R , CLC ,'CC L R ) 

390 CONTINUE 
Rfc TURN 


RfiGULAk INTERPOLATION BETWEEN TABLES 


ADO DO A LO K = NS»NP 

AID CCZ. 1(K)=TERP( C AM3K rCAMB ( K ) , C AMD T , CCZ 1 ( K ) ,CCZ2(K) ) 

RETURN 

END 

SUBRUU'i I NIL ARC (.RR AY , T AU , NROin S ,'NCOL S » I E t'WHERE » NLVLS ) 

SUBROUTINE TO INTERPOLATE BETWEEN LEVELS OF 
A GIVEN TABLE 

DIMENSION ARRAY (5, 25,12) ,NkGW S ( 6 ) , NCOL S < 6 ) » CLR ( 2 ) » ALR ( 2 ) ,TONY( 19) , 

1 WHERE (2, 6) »'YHERE(.2,'6) ,-MYCUL ( fc ) ,MAYY( 6) , CM ( 1 9 ) , XHERE < 2 , 6 > .MAWW (6 ) , 

2 MW COL (A) ,C(.19) ,’tPS(.19 ) , TRANS ( 1 9 ) , R F Y (. 1 9 ) , t T A ( 1 9 ) ,HOPP( 19) ,CLfAX< 19 
3) , /HERE I 2,6) ,HFRE ( 2,6) , Y( 19) ,UAU( 1 9 > ,'BET A ( 1 9 , 19 ) , TR I X < 19 , 19 > , MAZZ < 
Afa) , MZC0U.6) , [•; A X X (.6 )' , M X C 0 L (.6 ) 

COMKJN I N:\iUw, I S « IT (. 3 ) , A L PHA , RE. YN , C V A L , R t YUN , XMA X , ALMA X ( 19) ,CLMAX,C 

1 , EPS, TRANS, REY, 'ETA, HOPP , Z HERE ,'HE RE , ARRAY ,Y, BETA , TFAC,'TR1X,UAU,MAXX 

2 ,-MAZZ , MXCOL , MZCOL ,ASPEC , TAPER, 'BF , R F YND ,'D I SC R , P I ER , CRB , Q , T S T AX , EDGE 

3 * S IX , A.L PHR »'NF L AP , NL VL * N P , IY,1Z,IR,IP,'IV,ISTAR,A,B,H,TAUT,TAUR, 

A) w I-ST , R , BW X , YHERE ,-MYCUL , MAY Y , F L AP , T ONY »'TW1SA,X,Z. , CM , ACC , XHERE , 
5MWCUL , MAW« , CAMB ( 19) , CAM HR ,‘C AMBT , OUMY 1 , DUMY2 , NAME ( 25 ) ,AHERE<2,6) , 
6NAAA(6) , N /\C H L (•/ ) ,bhERt<2,6) » M A B B ( 6 ) ,MBC0L(6) ,CHERE(?,6) ,MACC(6) , 
7MCC0LI6) , CHERf (2,6) , M ADD ( 6 ) , MDC0L(6> 

I L = 0 

XTkAL=0. 

x r k a u = o . 

IS TAU LESS THAN UR EQUAL TO LOWEST LEVEL 
TAU VALUE 

1 B (. T A U— w H L R t: ( 2 , 1 ) ) 10, '20, 20 

LO J = 2 

X T R A L = 1 . 

GL1 TO UO 

20 IF (.T AU-WHEP.t- (2, NLVLS) ) AO, '30, 30 

YES,' SLT J ECUAL TO MAX LEVEL 

30 J = N L V L S 
XTkAO- 1. 

GO TO 80 

NU, SEARCH LIST QF LEVELS (TAU VALUES) 

UNTIL UNE GIVEN IS GREATER THAN ONE 10 BE 
FOUND 

AO l)0 50 J = 2, NLVLS 

If (.TAU-WHERE ( 2 , J ) ) 90, '90, 50 
50 CONIINUt 

ERROR IF NONE CAN BE FOUND 

I E = 6 

IF (.IE) 7G,70,t>0 

60 WR I r L ( I- P , 5 A 0 ) IE, 'INNOW, CVAL, ALPHA, RLYN,'REYON,XMAX 
70 RETURN 
80 CONTINUE 
90 L V L = w h t P. E ( 1 , J ) 

I" A X .< = . J (•. i.. W S (• J ) 

NAXC=NCOLS (J ) 
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I F CALPfiA— 999 • ).• 100fl30, ; 100 
100 I F ( REYN— 999 •! ) . 1 10 ,■ 120 ,• 1 10 
110 I F (XVA.L-999 • ) 260 ,250 ,<260 
120 1F6CVAL-999. ) 240F2 60 ,'2 4-0 
130 IPLCVAL-999. ) 140, 150, ‘140 

140 IF (REYN-999. ) 230, >260,^30 
150 IFUEYN-999. ) 260 ,■ 1 60 ,260 
160 IF (.R6Y0N-999. ) 1 8 0 1 70 ,• 1 80 

170 IF (.XMAX ): 2 2 0, 260 ,220 
180 IF UMAX) 200,190,200 
190 I F L P= 1 . 

00 TO 210 
200 I F L P = 2 
2.10 I S = 5 

GO TO 270 
220 1 S=4 

GO TO 270 
230 I S= 1 

GO TO 270 
240 I S = 2 

GO TO 270 
2 50 1 S = 3 

GO TO 270 
260 I F=6 

RETURN 

270 CALL LOOK ( ARRAY, L V.L , MAXR , MAXC ,IL ) 

CLR(2)=DUKY2 
AL R ( 2 ) =iJUMY 1 

GO TO (.280,290,300,310,320), IS 
280 F I ND = AL PH A 
ALPHA=999. 

GO 10 350 
290 F I ;nj D = K E Y iM 
R E Y N = 9 9 9 • 

GO m 350 
300 F I N U = C V A L 
C VAL=999. 

GO TO 350 
310 F I N 0 = R L Y Q N 
R E Y I J N = 9 9 9 * 

GO TO 350 
320 F I N D = X M A X 

GO TU (-330,340), IFLP 
330 X H A X = 0 • 

GO TO 350 
340 X M A X = 1 0 0 . 

350 LVL=aHlRE( 1, J-l) 

MAXR=NRUWS (.J-l ) 

M A X C = N COLS (.J-L ) 

SPECIAL PROCEDURE- FUR INTERPOLATION IN THL 
NEIGHBORHOOD OF CL MAX 

CALL LOOK ( ARRAY, L VL , MAXR , MAXC , It ) 

C LR ( 1 ) - o U M Y 2 
ALR ( 1 ) = D U M Y I 
IF ( IE) 360, *360, 60 
360 Kl=wHERL ( 2 ,-J-l ) 

R2=TAU 

K 3-= WHERE ( 2 J ) 

C3=F1ND 

GO TU (-4 9 0 ,*500, 370,620, 530)., I S 
370 I F ( INNUiv-1 ) 380,390,380 

380 I F (. I NNOK-5 ) 510,390,510 

390 OUT; Y1 = TLRP(.R1 , R 2 , K 3 , CL R ( 1 > , C L R ( 2 ) ) 

DUMY2=TlRP(.R3 ,R2,R3,ALR( 1) , ALR ( 2 ) ) 

1 f (■ ( ALPHA— AL R (■ 1 ) ) * L A L P 1 1 A - A L R ( 2 ) ) ) 400,5 1 0,5 10 
400 I F ( X T R AL - L.» ) 4 10,420,410 

410 If- (X TkALJ - L# ) 440,430,440 
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420 CVAL=OUKYL-CLR( 1 )+TfcRP (.ALPHA , UUM: Y2 , ALR (.1 > ,CLR( 1 ) ,CVAL) 

RETURN 

430 CVAL = DUMY1-CLR(.2 V. + TERP (.ALPHA , 0UMY2 , 'ALR (-2 >,CLR(2),C3) 

RETURN 

440 I F (• A L K ( 2 ) — A L R (. 1 ) ) 450 ,'450 ,‘460 

450 I F (.ALPHA— DUMY 2 ) 470,470,480 

460 I F (.ALPHA— DUMY2 ) 480,480,470 

470 CLA=TERP (ALR ( 1 ) , ALR ( 2 ) , ALPHA, CLR(1L»CVAL ) 

CLC = TERP (.R 1 , R2 , R3 ,'CL A ,'CLR L2 ) ) 

CVAL=TEKP( ALR (.2 ) , 'ALPHA , DUMY, 2 , CLC ,'DUMY 1 ) 

RETURN 

480 CLH = TERP(.ALK(2> ,ALR( 1 ), ALPHA, CLR(-2).,C3) 

C LC = T ERP ( R 1 , R2 , R3 ,'CL R ( L ) , CUB ) 

CVAL = TERP ( ALR (•! ) ,' ALPHA ,'DUMY2 , CLC , DUMY 1 ) 

RE TURN 

490 C 1 = ALPHA 

ALPHA=TERP(.R1 , R2 , R3 , C 1 , C 3 ) 

Rt TURN 
500 C1=REYN 

R G'YN= TERP ( R 1 , R2 , R3 , C i ,03 ) 

RETURN 
510 C l = CVA .L 

CVAL= I ERPLRl ,K2,R3,C1,C3> 

RETURN 

520 C 1 = RE YON 

RE'Y0N = TERP(.R1 ,'R2 ,'R3 , C 1 , C3 ) 

RE 1 URN 
530 C 1 = X M A X 

X H A X = T t k P ( R 1 , R 2 , R 3 , C 1 ,'C 3 ) 

RETURN 

540 FORMAT (-20X »■ 1 2HERRUR CODE (,I2,42H ) IE 1 CVAL GREATER THAN M AX VAL 
iUt LISTl U/37X,'34H1F 2 CVAL GKEATtR THAN TABLE VALUE/37X , 41H I F 3 AL 
2 P h A VALUE GREATER THAN TABLE V ALUL / 37X ,'29HI E 6 TAU VALUE CANNOT BE 
3EUU TU/3/X, I 31 i T ABLE NUMBER , I 2 ,■ 8 H IN CORE / 10 X , 8HC- V ALUE= , F I 0 . 2 , 5 X , 

4 6HA LPl i A= , E 1 U . 2 , 5 X »• L 31 IRE YNUL DS NO . = ,'F 10 . 2 , 5X , 6HRE YON= , F 1 0 . 2 , 5X , 5HXM 

5 A X = , E 1 0 . 2 / I 2 0 (• L H / ) ) 

Fill; 

SUBROUTINE LOOK. ( AR ,LVL , NAXR , MAXC , I E ) 


SUBROUTINE EUR TWO DIMENSIONAL TABLE LOOKUP 

D I MENS ION A( 5,2 5, 12 ) , TONY (. 1 9 ) , YHERE ( 2 , 6 ) , MYCUL ( 6 ) , MA Y Y ( 6 ) , CM ( 1 9 ) , 
1XHERE < 2,6) ,.\AWW(6 ) , MWCOL ( 6 ) ,C (■ 19 ) , EPS ( 19 ) , TRANS ( 19) ,REY( 19> ,ETA( 19 

2 ) , h C P P ( 10) ,CU*AX( 19) , ZHEREt.2,6) , HERE (.2,-6) ,Y( 19) ,UAU< 19 ) , BE T A ( 1 9 , 19 

3 ) , TRIX (.19, 19) ,MaZZ (.6) ,'MZCOL ( 6 ) , MAXX ( 6) ,-MXCOL ( 6 ) 

C 0 M w I J k I NNOw , I SB IT (.3 ) , AL PM A , REYN, CVAL , RE YON ,XMAX,ALMAX(19),CLMAX,C 
l, EPS, IRANS, 1 RE Y. , ETA, HOP P,ZT f ERE, HER L, 1 A, 'Y, 1 BETA, TFAC,TRIX,UAU,MAXX, 
2i\AZZ. ,MXCOL , KZCOL , ASP EC, TAPER, BE, RE YNL),D1 SCR, PIER, CRB, Ci , T STAX, EDGE, 
3SIG, ALPh R,NI LAP, NLVL ,NP, I Y ,■ I Z ,' I k ,' I P , I S ,• I STAR , A A , B , H , TAUT , TAUR , 

4 IBIS r,R, BEX, YE.ERE,MYCGL, MAY Y, FLAP, TUNY, Tw ISA, X,Z, CM, ACC, XHERE, 
BNwCOL , f AWu , GAMES ( 19 ) , CAMISk , CAMb T , DUMY 1 , DUMY2 

SECTION 10 LOCATE REYNOLDS NUMBER IN TABLE 

X T K A L = 0 . 

X TRaU = 0 . 

LROW IS ThE' LAST ROW CONTAINING INFORMATION 
IN ANY TABLE. THE LAST TWO ROWS OF ANY 
TABLE CONTAIN MAX DATA.’ 

L ROW = KAXR— 2 
IS=1 

IS RtYON EQUAL TO 999. 

IE YES THEN NORMAL AND WE WILL ALREADY HAVE 
SET THE SWITCH IS TO 1 AND WE -WILL USE THE 
VALUE OF RLY.N 
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IF NO T H t N SET REYN TO REYON AND SET SWITCH 
IS TC) 2 


If (REYUN-999. ). 10,20, 10 
10 REYN=REY0N 
I S = 2 

IS THE REYNOLDS NUMBER WE ARE LOOKING FOR 
GREATlK THAN OR EQUAL TO THE REYNOLDS 
NUMBER GIVEN ALONG THE FIRST ROW OF THE 
TABLE 

20 IF (i<EYN-A(LVL, l.MAXC) ) 40,30,30 
30 LOCR-MAXC 
X TR AU=1. 

GO 70 7 u 

NORMAL CASE.. IF REYNOLDS NUMBER IS LESS 
THAN MAX REY.N NUMBER IN TABLE THEN WE 
CHECK TO SEE IF IT IS GREATER THAN THE 
MINIMUM REYNOLDS NUMBER LISTED IN THE TABLE 

40 If (.RFYN-A ( L VL , 1,2) ) 50,60,'6C 

IF LESS THAN MINIMUM WE WANT TO USE THE 
MINIMUM VALUE IN THE TABLE 

SO L 0 C R = 3 
X T KAL= 1 . 

GO 10 70 

IF LESS THAN MAXIMUM AND GREATER THAN 

MINIMUM WE WANT TC SEARCH TABLE UinTIL WE 

FIND A REYNOLDS NUMBER GREATER THAN THE 

GIVFN REYNOLDS NUMBER 

IF REYN LESS THAN KEEP LOOKING 

IF REYN GREATER THAN SET UP R1,R3 FOR 

INTERPOLATION 

60 DO 610 LGCR = 3 » MAXC 

IF (.KEYN-AI LVL ,'l , LOCK ) ) 70,70,610 
70 R 1 = A ( LVL , L ,'LOCR - 1 ) 

R 3-= A ( LVL , 1 , LOCK ) 

GO TO (HO, ISO), IS 

SECTION F UK NORMAL LOOK UP OF ALPHA AND 
OVAL. A VALUE OF 999 SPECIFIES THAT THIS 
VARIABLE IS THE ONE whUSE VALUE IS TO BE 
FOUND IN THE TABLE. 

80 IF ( ALPHA-999.. 1 90,110,90 
90 IF (.KEYN-999. ) 100, 600, '100 

100 IF (.CVAL-999. ) 600 ,190 ,'600 
110 I F ( C V A L - 9 9 9 . ) 120, 130,120 

120 1 F (.REYN-999 . ) 160 ,'600,160 

130 IF (.REYN-999. ) 600, 140, >600 
140 IF ( KCYUN-999 . ) 150,600,150 

SECTION FOR LOOKING UP MAX VALUES 

150 IF UMAX) 580,590,580 

SECTION TO LOOK UP ALPHA FOR GIVEN CVAL BY 
INTERPOLATING BETWEEN COLUMNS WHERE 
REYNOLDS NUMBERS BRACKET GIVEN REYN 

L60 on 170 J = 3 » L R 0 W 

CCC-TERP (-K 1 , RE YN ,'R3 , A ( 1 VL J , LOCR-1 > , A ( LVL , J ,LOCR ) ) 

IF (C V A L-CCC ) 180,180,170 
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170 CONTI NUt 
1 F = 2 

GO T (J 600 

180 CC2= IL-RP (R1 * KE YN , k 3 , A ( LVL ,'J- 1 , LOCR- l >. , A ( L VL , J- 1 , LOCR ) ) 

ALPHA = T C RP LCC2 , C VAL » CCC * A ( L VL » J— 1 , 1) , A <.L VL,-J , 1 ) ) 

CO TU 600 

LOOK UP CVAL FOR GIVEN ALPHA. CHECK. TO SEE 
IF WE HAVE A LIFT OR DRAG TABLE IN AT THIS 
TINE. DRAG TABLES HAVE ZERO VALUES FUR ALL 
MAXI-NUN VALUES (.LAST TWO ROWS). 

LOO IF (.A ( LVL ,MAXR,'LOCR-.l ) ) 210,-200,210 

200 IF (. A ( L V L , M A X R »' L I ) C R )! ) 210,550,210 

FOR LIFT TABLE. DUMY1 AND DUMY2 ARE THE 
MAXIMUM INTERPOLATED VALUES FOR GIVEN 
REYNOLDS NUMBERS OF ALPHA MAX AND CVAL MAX 

210 CONTINUE 

DUMY1 = T ERROR 1, REYN, 'R3, 'A ( LVL, MAXR,L0CR-1 ) , A <-L VL , MAXR , LOCR ) ) 

DUMY2 = TERP (.R 1 , RE YN 3 ,'A ( LVL , MAXR- 1 ,'LOCR- 1 ) , A ( LVL , MAXR- 1 , LOCR) ) 
1F(XTRAL-L. ) 220' ,'2 A 0 ,'220 
220 IFUTRAU— L. ) 230 ,'2 AO ,'230 
230 CONTINUE 

I ft ( AL PI tA-A< LVL, MAXR, LOCR- 1) )*{ ALPHA-A (.L VL , 'MAXR , LOCR ) ) ) 2 AO , 2 50 , 250 
2 AO IF (.A (LVL, MAXR, LUCK) -A (LVL, MAXR, LOCR- 1) ). 270,260,260 
250 GO TO 560 

260 ALPH 1 = A ( LVL ,MAXk, 'LOCK- L) 

ALPIl2 = A( LVL , MAXR , LOCR) 

GO TO 20 0 

270 A L P h 1 = A ( LVL, MAXR ,-LOCR ) 

ALPH? = A( LVL ,"MAXA, LOCR-L ) 

280 IF (.XT RAL-L. ) 290 ,'300 ,'290 
290 IF (.XT RAO- L. ) A GO ,'390 , A80 
300 I F (.ALPF'2— OUMYL ) 350,550,310 

310 IF (.( ALPF.A-DUMYl ) * ( ALPHA-ALPH2 ) ) 320,550,550 
320 DO 3 30 J - 3 » LRCw 

IT (.ALPH2— A ( LVL, J, 1) ) 3A0,3A0,330 

330 CONTINUL 

I E = 3 

GO TO 600 

3A0 C 1- A ( LVL,J,'L(ICR-l ) 

CVAL = Dl.P Y2 + ( A (.LVL ,'MAXR- 1 , LOCK-1 >-C 1 ) / ( AL PH1-ALPH2 ) * ( ALPHA-DUMY1 ) 

GO TU 600 

350 IF ( (ALPHA— DUf-Yl )*(ALPHA-ALPH1) ) 360,550,550 

360 uO 370 J=3 ,LkOW 

I F ( ALPI. 1— A (.LVL ,-J , 1 ). ) 380,380,3 70 

3 70 CONTINUL' 

IE =3 

GO (U 600 

380 C1 = A( LVL , J , LOCR-1 ) 

CVAL=UUMY?+ ( A (.LVL ,'MAXR- 1 ,L0CK-1 ) -C 1 ) / ( AL PH2-ALPH 1 ) * ( ALPHA-DUMY 1 ) 

GO TO 600 

390 IF (.ALPFI2-JUMY1 ) A00,5S0,AA0 

AOO I F (. ( ALPHA— DUMY 1) M-ALPHA-ALPH1 ) ) AID, 550, 550 
AID DO A 2 9 J = 3,'LR0W 

If (.ALPH 1 — A (.LVL, J , 1 )) A30, A30, A 20 
A20 CONTIimUE 
I L= 3 

GO TO 600 

A 30 C 1 = A( LVL , J , LOCR > 

CVAL-=l;Ul v, Y2+ ( A (.LVL ,'MAXR— 1 , LOCR >-C 1 ) / ( ALPH2-ALPH1 ) » ( ALPHA-DUMY1 ) 

GO TU 600 

AAO IF(. (ALPUA-DUf-Yl )*( AL PH.A- ALPH2 ) ) A50,550,550 
A 50 on A60 J = 3 ,' L R 0 W 

IF (.ALPH2-tA(.LVL,J,'1 ). ) A70, A70,A60 

A 6 0 CON T I NUF 
I L = 3 


161 


oonomonnrooono ooo 


GO ro 600 

470 Cl=A(jLVL,J,L0CR) 

CVAL = DUMY2+ ( ALLVL,'MAXR-l,LCCR)-Cl )/( ALPH2-ALPH1 )*( ALPHA-DUMY1 ) 
GO TO 600 

480 I F (ALPHA— 0UMY 1 ) 490,490,5 20 
4)90 DO 500 J-3 ,'LROiV 

IF (.ALPH1— A(.LVL, J, 1) ) 5.10, 5L0, 500 
500 CONTINUE 
I 6 = 3 

GO TO 600 

5 LO C l=A (1 VL , J , LCiCR-1 I 
C3=A(LVL,J,'L0CR) 

Cl=TERP!Rl,REYN,R3,Ci,C3) 

CVAL=TERP ( A L PH 1 , ALPHA ,<OUKY L , C 1 r 0UMY2 ). 

GO TO 600 

520 DO 530 J = 3 ,-LRGW 

IF LALPH2— A (.LVL, J , 1 ) ) 540,540,530 
530 CONTINUE 
I E= 3 

GO TU 600 

540 C1=A(LVL,J»L OCR-1) 

C3=A(LVL*J,L0CR) 

C3=TEkP(R1,'REYN,R3,C1,C3) 

CVAL = TERP ( DUF Y 1 , ALPHA, ALPH2 , DUMY2 , C3 ). 

GO TO 600 

550 DO 560 J=3 , LkOW 

SLARCH FOR ALPHA 

IF (.ALPHA— A (.LVL , J , 1 ). ) 570,5 70,560 
560 CONTINUE 
I E= 3 

GO !U 600 

570 CCC = TERP(Rl,REYN,R3,A(LVL,J-l,L0CR-il ) , A ( LVL , J-l , LOCK ) ) 

CC2 = TtRP (.Rl ,KtYN,R3, A( LVL, J , LOCR-1 )., A ( LVL , J , LOCR ) ) 

CVAL = TERP(.A( LVL, J-l, 1 )., ALPHA , A ( LVL , J ,'l ).,CCC,CC2) 

GO TO 600 

580 XMAX=TERP(R1 , REYN, i<3, ALL VL,KAXR, LOCK- 1 ) , A ( L VL , M A XR , LOC R ) ) 

GO 10 (600, '620), IS 

590 XMAX=TERP(.:U-, -REYN, R3, -A (-LVL, MAXR-1, LUCK-1 ) , A ( L VL , KAXR- 1 , LOCR ) ) 
GO TO (£.00, '620), IS 
600 RETURN 
610 CUNTINUt 

1 E= 1 

GU TO 600 
620 REYN=999. 

GO TlJ 600 
END 

OVERLAY ( 6L INDA, 1,0) 

PROGRAM ON E 
CALL RAIN 
RE TURIN 
E NO 

SUBROUTINE MAIN 


SUBROUTINE MAIN IS THE MAIN PROGRAM FOR 
AERODYNAMIC ANALYSIS OF STRAIGHT WING 
AIRCRAFi. THIS IS THE FIRST OF A SERIES OF 
OF MAIN SUBROUTINES. 

THE SEQUENCE IS 

MAIN 
M A I N I 

(FLAP) IN ON FLAP) 

NAIN3 MAIN2 


N A I IN' 5 


MAIN4 

LEND OF PROGRAM) 


L 1 F f 1 WHERE MXCUL MAXX 
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DRAG 1 

ZHERE 

MZCOL 

MAZZ 

PM 1 

YHERE 

MYCOL 

MAYY 

PMM1.F LAP ) 

XHERE 

MwCOL 

MAWW 

LIFT 2 

AHERE 

MACOL 

MAAA 

DRAG2 

BHERE 

MBCOL 

MABB 

PM2 

CHERE 

MCCOL 

MACC 

PMM2CFLAP) 

DHERE 

MDCOL 

MADD 


DIMENSION C ( 19), EPS ( 19 ) , TRANS (.19 ) ,KEY( 19 ) , £ T A ( 19 ) , HOPP ( 1 9 ) ,CLMAX< 1 
19 > t'ZHERE (.2,6 ) .WHERE ( 2 ,'6 ) , ARRAY ( 5 ,'25 , 12 ) , Y ( 1 9 ) , TAU ( 19 ) , BETA ( 19 , 19 ) , 
2MAZZ ( 6 ) , MZCOL (.6 ) ,'MAXX(.6) , MXCUL ( 6 ) , YHERE ( 2 , 6 ) ,MYC0L(6) , MA YY ( 6 ) ,TRIX 
3( 19,19) , CM (.19):, DISK 19):, YD A ( 19) , YDAX ( 19) ,YX( 19) , TONY ( 19 ) , ALPG ( 1 9 ) , 
A A U P H Z ( 19) , X HERE ( 2 ,-6 ) ,'MWC0L(.6) ,'MAWW(.6)! 

COMMON INNOW, ISWITI3) , AL PHA , RE YN ,'CL L ,'RE YON ,'XMAX , ALMAX < 19) ,CLMAX,C, 
1EPS, T1;ANS, KEY, iETA »«HOPP,'Z HERE, WHERE, 1 ARRAY, Y,tBETA,TF AC, TR IX, TAUtMAXX 

2 , MA ZZ, MXCOL, MZCOL, ASPEC, TAPER, BF, 'RE YND.'D I SCR, PIER, CRB, Q.TSTAX, EDGE 

3 , S I-G , AL PHR ,NF LA P ,'NL VL ,'NP , I Y. , I Z , IR, IP, IS, I ST AR , A , B , H , TAUT , TAUR , 
ATWI-ST, R,BWX,YHERE„MYCnL,MAYY, FLAP, TON Y ,'T W I S A , X , Z , CM , ACC, XHEKE, 
5MWCOL , MAWW , CAMB ( 19 ): , CAMBR ,'CAMBT , DUMY 1 , DUMY2 .NAME ( 25 ) ,'AHERE (2,6), 

6 M A A A ( 6 ) , MACOL (.6) ,BHERE(.2,6 )., M.ABB ( 6 ) , MBCOL ( 6 ) ,CHERE( 2,6) ,MACC<6) , 
7MCCOL ( 6 ) ,UHbRE(.2,'6) , MADD ( 6 )! , MDCOL ( 6 ) 

READFR (IR) AND PRINTER (IP) LOGICAL UNIT 
NUMBERS 

I R = 8 
I P = 6 

RF„ [NO 1 
KE W I NO ? 

REWIND 1 
REWIND A 
R T V) I N D 6 
REWIND 1C 
REWIND 15 
REN I In LJ 20 
REWIND 7 

INPUT DATA SECTION 

IS=1 

REALM IR, 6 60) ASPEC ,'T AUT , TAUR , TAPER ,'TW I ST ,R,.BF , RE YND ,DI SCR , A , B , H , 
1ALPHK,N(-LAP,HAP,'X,7,'TWISA,'CAMBT,'CAMBR 
I F (AS PEC-99 . ) 20,10,10 
10 CALL EXIT 
20 CONTI NUL 


LAYOUT 

OF 

FOURTH DATA CARD 


FIELD 

1 

I 1 

NUMBER OF TAU 

VALUES PER TABLE 

FIELD 

2 

I 1 

1 FOR READ IN 
0 FOR NO READ 

OF TAU, RE Y , E TC . 
IN 

FIELD 

3 

I 1 

1 FOR DUMP OF 
0 FOR NO DUMP 

COMPUTED ARRAYS 

FIELD 

4 

I 1 

1 FOR DUMP OF 
0 FOR NO DUMP 

BETA ARRAY 

FIELD 

5 

ri 

1 READ CUBE 1 FROM TAPE, LOAD 
Tl) DISK, COPY CUBE 1 TO CUBE 


2 ON DISK. 

2 READ CUBE 1 FROM TAPE, LOAD 
TO DISK, READ CUBE 2 FROM 
TAPE, 'LOAD TO DISK. 

3 READ CUBE 1 FROM CARDS, LOAD 
TO TAPE, LOAD TO DISK, COPY 
CUBE 1 TO CUBE 2 ON DISK 

A READ CUBE 1 FROM CARDS, LOAD 
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TO TAPE i LOAD TO DISK, READ 
CUBE 2 FROM CARDS, LOAD TO 
TAPE,* LOAD TO DISK 

FIELD 6 25A2 FIFTY COLUMNS OF I DEN T I F I MG 

INFORMATION. THIS IS PRINTED 
AT THE TOP OF EACH PAGE OF 
OUTPUT 

READ (IK, 660) NLVl , I SW IT , IG f NAME 

• SWITCH ZERO ON WILL READ IN VALUES FROM 
CARDS FOR T A U REY, C, EPS, EDGE, CRB, AND 
ACC. FORMAT IS 16F5..0. ARRAYS ARE READ IN 
ROWWISE 

CALL SETSW 
DO 30 J = 1 , 6 
WHERE ( L , J ) - 0 . 

WHERt:(2,J)=0. 
to A W W ( J ). = 0 
30 M A C U L ( J ) = 0 
I Y - 0 

ITMG-3) 80, 40, AO 
40 DO 60 J=l,4 

CALL AERDA(.ARRAY,NLVL*yiHLRE,MAC0L,MAWW»7 , IY ) 

50 WRI TE ( 7 ) W HERE, MAC OL, RAW W 
IFUG-4) 80,60,80 
60 DO 70 J = 1 , 4 

CALL AERDALAKRAY, NLVL, 'WHERE, MACOL ,MAWW, 7, IY ) 

70 WRITE ( 7 ) WHERE, MCOL ,'MAWk 
80 I K = I 

REWIND 7 

REAUm ARRAY 

READt 7 ) WHERE ,'MXCGLt MAXX 

WRITE (IK). ARRAY 

I K = I K + L 

K E A D ( 7 ) ARRAY 

READ (7 ) ZHE RE ,'MZCOL , MA.Z Z 

WRITE ( IK) ARRAY 

I K = I K+ 1 

R E A D ( 7 ) ARRAY 

READ ( 7 ) YHERE ,'MYCCL , MAY Y 

WRITE(IK) ARRAY 

I K = I K+ 1 

READt 7 ) ARRAY 

READ (7 ) X H E R E M W C 0 L , M A'W W 

WRITL(IK) ARRAY 

REWIND I 

REWIND 2 

REWIND 3 

REWIND 4 

I F {- 1 G — 1 ) LOO, '100,-90 
90 I F ( I G — 3 ) 130, ICO , 1 3u 

LOO I K— I 

DO 1 L 0 J = L v4 
READt LK) ARRAY 
KK = IK* 5 
I K= I K +.1 

110 WRITE(KK) ARRAY 
DO 120 J = 1 , 6 
AHEKEU, J)=WHERFL1, J) 

A H E R E ( 2 , J ) - W H E R t (. 2 , J ) 

BHE RE ( L , J ) = Z HERE (. I , J ) 
bHEkE(2, J) = ZHERE (-2, J ) 

CHERE( L, J)*YHEKL-( 1 , J ) 

CHERE ( 2 , J ) - YHERb ( 2,J ) 

DHEKE ( L, J ) - X H E R b (• I , J ) 

□ HERE ( 2 , J ) = X H E R E (. 2 , J ) 
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f'lACUL ( J ) =MXC(JL C-J ) 

MB COL < J ) = M Z f, 0 L ( J J 
MCCL-L < J ) -MYCUL (.J ) 

MDCOU J)=MWCGL!J) 

MAAAI J ) = M A X X ( J ) 

M A 6 B ( J ) = M A Z Z ( J ). 

MACC! J)=KAYY{ J) 

120 foAOO ( J ).= MAWW ( J ) 

GO VC 1 A 0 
130 I K = 5 

HE AD < f ) ARRAY 
READ { 7 ) AHERL ,MACOL»MAAA 
WRITE! I-K > ARRAY 
IK= IK + 5 
RE AD ( 7 ) ARRAY 
KEaO( 7 ) BHERE , MBCU.L , N ABB 
WRITE! IK) ARRAY 
IK- IK + 5 
READ! 7 ) ARRAY 
RE AO ! 7 ) CEiERfc r MCCOL , MA.CC 
WRITE (I-K) ARRAY 
IK=IK + ‘j 
READ! 7) ARRAY 
READ! 7 ) DPERt, MCCOL, MAOD 
WRITE! IK) ARRAY 
140 REWIND 5 
REWIND 10 
REWIND 15 
R E W I iv D 20 
REWilNL 1 
REWIND 2 
REWIND 3 
REWIND 4 
ALPS! 1 ) =0 . 

NP=R-1. 

J P = ■'< / 2 . 

P I E R = 3 . 1A 159/k 
CALL DA T S w (-0 » I) 

GU TO (.15 0, '160) 


150 READ! I R , 6 70 ) 
READ! Ik, 670) 
RCA D < 1 R , fc 70 ) 
READ! I R , 6 70 ) 
READ! IR, 6 70) 
READ! I R , A 70 ) 
READ! I R , t>70 ) 

160 cocrioot- 


I 


( TAU! I ) , I = 1,NP) 

(.KEY!, I ) , I = 1,MP) 

!.C( I) ,* I = 1 » hi P ). 

(EPSU) , 1 = 1, NP) 

EDGE 

CRU 

ACC 


II NO FUSE (TUSELAGE) 

IE (A) A 10, A 10, 170 
170 YO=b*SCjKT(. l.-H**2/A*»2) 

ECC = SGRT <.A**2-T,**2 ) 

CALI DA V Sw (.0 , I ) 

GO TD ( 190,160) , I 

180 CRB = 2 . * < l.-YC)/ ( ASPEC*(.1.*TAPER-2.*Y0*TAPFR) ) 
190 TFAC=l.-(.Y0*TAUk*CRB/(.3.1ALB9*A*B) ).*4 . 
j P j P ^ j 


Y PRIME ( I ) 

DO 200 I=1,JPP 
XI I=i-L 

200 Y.DA! I )=Y0+!.1.-Y0)»C0S!XI I*PIER) 

CHECK IF FUSE ELLIPTICAL OR CIRCULAR 

IF(A-D) 230,230,210 
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ELLIPTICAL FUSE 


210 CONTINUE 

DI STX=0.5*(.SQRT( 1 . + ( H-CCC ) **2 ) +SQRT ( 1 . +. ( H + ECC ) * * 2 ) ) 

BWX=l./( A— 8)*(.A-ti*.DISTX/SQRTU)ISTX«;*2-ECC**2) ) 

00 220 I=*.1,JPP 

U 1 ST ( I ) =0. 5* ( SORT (. y 0 AU )*» 2 +( H-ECC ).**2 ) + SQRT ( YD A ( I ) **2+ ( H + ECC ) **2 
1 ) ) 

Y BAR PRIME (I) 

220 YD AX ( I ) = UDA( I ) / l.A-B ) * C.A-ti*.D I S T ( I ) /SORT (DISK I )«*2-QCC*»2) ) ) /bWX 
CO TU 2 fj 0 

Y BAR PRIME (I) FOR CIRCULAR FUSE 

230 DO 24 0. l = l,'JPP 

240 YDAX ( I ) = YCA ( I ).* ( I .-A* *2/ ( YD A ( I ) **2+.h**2 ) ) 

B WX = 1 . - A * # 2 / ( 1 ♦ + H * *. 2 ) 

COMMON TO ELLIPTIC AND CIRCULAR FUSE 

250 DO 260 I=!L,JPP 

Y BAR (-1 ) 

A I = I - I 

260 Y X (• I ) = C 0 S (• A I * P l E R ) 

Y (.1) 

DO 270 1=2, JPP 

270 Y (. I ) = 1 EkP ( YDAX (. I - 1 ) , Y X (. I ) ,'YCAX ( I ) , Y.DA ( I - I ) , YDA ( I ) ) 

DO 200 I = 1 » J P 
Y (. I- ) = Y (• I +. 1 ) 

280 Y X (• 1 ) = Y X (. 1 + i ) 

M = J P - 1 

DO 290 1 = I M 

IR1=R 

I R i= I R 1 - l 

YX (. I R I )=-7X( 1 ) 

290 Y.MRI )=-Y( I ) 

If(A-h) 360,360,300 

FUSE IS ELLIPTICAL 

300 00 310 I = L ,■ J P 

DIST ( I )=0.5* ( SORT ( Y( I ) .*.*?+. (.H-bCC ) **2 ) +SQRT (Y ( I ) * *2+ (H + ECC ) **2 ) ) 
Q=U 1 ST ( I )'/SORT(.i..;IST( I )**2-ECC**2) 

S=l. + ( ECC*Y ( 1 ) / ( DISK I) **2-LCC*»2 ) )>*2 

1 R A 'J S ( I ) = L./ ( A-li ) * (. A - B *. 0 / S ) 

I R 1=K 

I R I = I K I- I 

310 1KANSI I R l ) = TRANS (• I ) 

WING UN TOP OR BOTTOM OF FUSELAGE 

IF(A-H) 320,320,330 
320 TRANS <JP)=1. 

330 IF(.NELAP) 5 20 ,320 , 340 
3A0 IFUE-L.) 350,520,520 

350 DI Sl=0. 5*( SORT {.t;F**2+ ( EI-LCC ) **.2 ) +SORT (BE**2+(H+F.CC)**2) ) 

BFX=( tiF / (.A-h ) * {.A — 13*1) I SE/SOR F( D I SF**.2-ECC **2 ) ) ) /BwX 
GO TO A 60 
360 DO 370 I — L , J P 

Tl<ANS( I ) = 1 . +A**.?« (Y ( I )«*2-H**2)/( ( V ( 1 ) ** 2+H* * 2 ) * *2 ) 
i r t=k 

I R I = I R I - I 

370 TRANST IRI ) = fKANS(.I ) 
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IF(.A-H) 360 * 3 80 f 39 0 
380 TRANS ( JP ) = 1 • 0 
390 1 F (.61 — 1. ) 400,520,520 
400 BFX=( BP* ( I • - A*.*,2 / (,BF * *2 +H* *,2 )).)/BWX 
GO TO 460 
4.10 00 420 1 = 1, NP 
A 1= 1 

Y {.I- ) =C 0 S (. A I * P I E R ) 

YX (. I ) -Y ( I ) 

42.0 TRANSU )=*1. • ■ . 

YO— 0 •' 

SPECIAL CASE IF VALUES HAVE BEEN READ IN 
DO NOT WANT TO COMPUTE CRB • 

WILL NOT COMPUTE VALUES ALREADY READ IN 

CALL DATS W TO , I ) 

GO TO (-440,430), I 

430 CRB = 2 . * ( l.-YU)./(ASPEC*U.-hTAPER-2.*YU*TAPBR> ) 

440 T F AC = 1 . 

B.W X = 1 • 

IF6BE-1. T 450,520,520 
450 BFX-BF 

460 TSTAX = ACOS(-bFX/BWX ). 

DO 510 1=1, JP 
A I = I 

C3=AI#PIER 
TSTX=TST AX 

IF(.C3 — TSTX) 510,470,460 
470 FSTAX=C3 
I S T A K = I 
GO TO 520 
480 A M = C 3 - T STAX 

Cl=( Al-1.. ) * P I E R 
AK=TSTAX-C 1 
IF (. AM- AK ) 490 ,490 ,500 
490 ISrAR=I 
T S T A X = C 3 
GU TO 520 
500 ISPAR=I-1 
rSTAX=Cl 
GO TO 520 
510 CONTINUE 
520 CALL DATSWLO, I ) 

GO TU ( 5 5 U , 5 3 0 ) , 1 

5 30 EDGE = SORT 6 1 . +4 . / ( A S PLC* * 2 ) ) 

DO 540 1 = 1, *!JP 

C (. I )=!.-< L.- TAPER )MABS(Y(.I ) ) -YG ) / (.1 . -YO ) 

TAU ( 1 ) =T AUK/C (.1) * (.1 .- ( 1 1 APLK# TAU1 /TAUk ) * (.ABS ( Y ( I ) ) - YU ) / ( 1 .-YO) ) 
CAM B ( 1 ) = CAMi*R+.(.CAMrtT-CAMBK )'* ( AB S ( Y (. I ) ) - YC ) / ( 1 •- YO ) 

ACC = 0.666 6* ( 1 . +1 APfcK* ( 1 . + TAPt K ) ) /( 1 . + TAPER). 

540 KFY( I ) = R E Y N D * C (■ I )/ACC 
550 L0CER=3 

A L P H A = 9 9 9. 

C LL = 0 * 

R LYON =999. 

X M A X = 0 • 

LOOK UP ZERO LIFT ANGLE FOR TIP SECTION 

CALL t'K I DG (. A L PH7 ,1,1, 1 ,1 Y, 1 IL , RE Y , TAU»-CAMB , ALPHA , 1 , 1 , ALPHA) 

1 T ( 1 1 ) 560 ,'560 , 640 

LOOK UP ZERO LIFT ANGLE FOR ROUT SECTION 

56C A L P El A = 9 9 9 . 

CALL OR I OG (« A L PUZ ,■ J P , J P , 1 , I Y , 1 E , RE Y ,-TAU ,-CAKB , ALPHA , 1 , 1 , ALPHA ) 

I F { Eh) 570,-57 0 ,64 0 
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CALCULATE AERODYNAMIC AND GEOMETRIC TWIST 


C 

c 

c 

c 


c 


c 

c 

c 


c 

c 

c 


570 IT ( TWISA-1C0. ) 580 ,590 ,590 
580 TWIST = TwISA-(ALPHZUP)-ALPHZ< 1) ) 

GO Tfj 600 

590 TW I-SA= TWIST+(ALPHZ(.JP) — ALPH7 ( i) ) 

600 CALL OATSW (.0 , I- ). 

GO TO (550,610), I 
6.10 DO 620 1 - 1 ,NP 

b20 EPS ( I ) = I Wl-ST*TAPEIs*( AtiS ( Y (. I ) ) -Y0 ) / (.C (• I ) * ( l.-YO) ) 

650 LOC ER= 1 
CLL=999. 

ALP HA =999. 

KLYN=999. 

RF Y0N=0 . 

FOR XMAX-100. ARC WILL LOOK UP ALPHA MAX 
FOR X M A X = 0 . ARC WILL LOOK UP CL MAX 

X M A X = 1 0 ( . 

CALL DR IDG CALM AX, 1 , r JP , 1 , I Y , I E , K L Y , T AU , C AMB ,AL PH A , 1 , 1 , XMAX) 

C A.L L M A I N 1 
RETURN 

640 WRITE ( IP, 680) IF.LOCER 
CALL EXIT 

660 FORMAT ( HE 10 . 0/.5F 10 . 0 , I 1 0 » 4F 5 . 0 / 2 F 10.0) 

660 FORMAT (.61.1, 25A2) 

670 FORMAT (16F6 .0 ) 

680 FORMAT MX, 1 IHERROr CODE ,I2,13H -IN SECTION ,13, OH ABORTED) 

END 

SUBROUTINE MAIN1 


CONTINUATION OF SUBROUTINE MAIN 

DIMENSION C ( 19 ) , EPS ( 19 ) , TRANS (.19) , RE Y ( 19 ) , E T A ( 1 9 ) , HOPP ( 19) ,CLMAX( 
119), ZHEkE (2,6).,WHLKt( 2,6 ) ,ARRAY ( 5, 25,12 ) , Y( 19 ) , TAU( 19 ) , BETA ( 19 , 19) 

2 , M A Z Z ( 6 ) , Mi Z C 0.L (.6 ) ,-MAXXU ) ,-MXCOL ( 6 ) ,-YHERE ( 2,<6 ) ,MYC0L(6) , MA Y Y ( 6 ) , 
3TRI-X (.19, 19 ) ,THE ( 19)., TONY ( 19) , GENE ( 19). , SIGMA ( 1 9 ) , LL ( 19 ) , MM ( 1 9 ) , CM ( 

A 19) ,XHERE(2,fc):,MWC0L(6) ,MAWW(6) 

COMMON INNOW, ISWIT(,3> , AL PHA , RF Yn , C LL , RE YON,'XMAX , ALM AX ( 19) ,CLMAX,C, 
1EPS, TRANS, RE Y, 'ETA, HOPP, ZHERE, WHERE, 'ARRAY, Y, 'BETA, TF AC, TR IX, TAU.MAXX 
2, MAZZ,MXCOL,MZCUL,ASPEC, TAPER, 'BE, REYNOLD I SCR, PIER, CRB, Q, TSTAX, EDGE 

3 *'S I-G, AlPHk , NFLAP ,NLVL ,-NP , I Y , I Z , I R ,• I P ,■ I S , I STAR , A , B , H , TAUT , TAUR, 
ATWIST,R,BWX,YHERE,MYCOL,MAYY,FLAP,TONY,TWlSA,X,Z,CM,ACC,XHERE, 

5 MW COL , MAWw ,'C AM.B ( 1 9 ) , C AMBR , C AMBT , DUMY 1 , DUM Y2 , NAME ( 25 ) , AHERE (2,6), 

6 M A A A (6 ) , MACOL (-6 ) , BHERE (.2 , 6 ) , MABB ( 6 ) , MBCOL ( 6 ) ,CHERE(2,6) , MACC(6) , 

1 MCCOL ( 6 ) , L'HhRE (2 ,6 ) , Mi A D D ( 6 ) , M.DCOL ( 6 ) 

COMPUTATION OF BETA 


DO 80 M= I , NP 
DO 80 K = 1 , N P 
AK=K 
A M = M 

IF(K-Mi) 20,10,20 

10 B E T A (• M, , K ) = 1 8 0 . * R / (. 8 . * 3 • 1 4 l 59 *S 1 N (AK*P I ER ) ) 

GO TO 50 

20 KK=IABS(K-M)-2«( (. I ABS (..K-M ) ) / 2 ) 

IF(-KK-L) 40, 30, AO 

30 BETA (M,K ) = lbO./ (4. *3. 14159*R*S IN ( AK»PIER ) >* ( 1./ ( l.-COS ( ( AK + AM) * 
13.14169/K). )-l./( l.-COS ( (AK-AM). *3.14159/R) ) Y. 

GO TO 50 
40 BETA (• E. , K ) =0 . 

50 IF(M-K) 70,60,70 

60 TRI-X ( i*1 , K ) = 1 . +.0 • 1*C (.K ) *CRB/BWX*.B L T A ( M , K )'* 1 RANS ( K ) 

GO TO 80 

70 TR I X ( M » K ) =0 • 1 *.C ( K ) *CRB/.BWX*BE TA ( M , K )* TRANS <-K ) 

80 CALL DAT Sw (.4,1) 
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GO TU (90,100), I 
90 WRITE ( IP, A 30) 

CALL SSSUlETA.NP) 

loo on no i=u,np 

A I — I 

110 H T A ( 1 )=(;.523598/R*(.3.-(.-l. ).**I ) * S I M ( A I * P I ER ) 

CHECK IF FLAP CASF 

IF(.NFLAP). 3 10, *3 1C, -120 

IS THERE A PART-SPAN TLAP 

120 IF ( OF— 1 • ). 130, '310, 310 
130 00 170 I - L ,' N P 
A 1 = I 

THE ( I ) =A I *P I cR 

IF ( THE ( I )-l STAX ) 150,140,150 

1 A 0 T0NY( 1 )=i./45. #TS1 AX*S 1N< THE ( 1 ) ) 

GO TO L6C 

150 TUNY( I-) = l./90.*( (-COS( THE ( I )!)-C0S( TSTAX ).) *(ALOG( 1 .-COS (THE ( 1 )+TSTAX 
1 ) )-ALOC( 1 .-CHS ( THE (. I )-T STAX ) ) ) +2 . *TSTAX* S IN ( THF ( I ) ) ) 

160 G E N F ( I ) = 1 . / 9 0 . * ( (-COS ( THE ( I ). ) +COS ( TS TAX T ) * ( ALOG ( l.+COS(THE< I ) -TSTAX 
1) ) -ALOG ( 1 . +.CCS (. I HE (■ I ) +TSTAX ) ) ) +2 . * (.3 .' 14 1 59-.T STAX ) *S I N ( THE ( I ) ) ) 

170 CONTINUE 

CHECK FOR UUN-P 

CALL DATSW (.3,1) 

GO TU ( 180,190) , I 
180 W R I- T E ( IP, -440) 

CALL AAA (. GENE, NP) 

WRI re < I P , 450 ) 

CALL A A., (.TONY , NP ) 

190 00 230 K - 1 ,• N P 
SUN’ = C . 

00 200 H- 1 ,'NP 

200 SUN =TONY (H ) *l'.E T A (>. , K ) +.SUK 

1 F ( K - IS r AH )! 210,210,220 
210 HOPPIK ) = 1 .-SUN. 

GO TO 230 
220 HOP P ( K ) =— S U 
230 CONTINUE 

DO 270 K - L * NP 
SUN = 0. 

DO 240 H - 1 ,’ N P 

240 SUP = GCNC (.M ) *BETA(.h, K ) + SUM 
I R I = K 

IF (.K-MI + ISTAR) 250,260,260 
250 S I G N A ( K ) - 1 . - S 0 N 
G(J II: 2 70 
260 S I G h A ( K ) = -SUP- 
270 CONTI N(JE 

Ul) 2b0 K = 1 ,'NP 

HOPP ( K ) = S IGN A (.K ) - HOPP ( K ) 

280 TONYIK ) =GENE (K)-TONY(K) 

CHECK FUR DUMP 

CALL UAT Sw (.3 , 1 ) 

GO TO (290, '300), I 
290 HR 1-TE ( I P ,'A60 ) 

CALL AAACHOPP.NP) 

WK I Ttf ( IP, 470 ) 

CALL A A A (. S 1 G M A , N P ) 

HR ITU IP, 440 ) 

CALL AAA(GFNE,NP) 

WRITE ( IP ,'450 ) 

CALL AAA(.TONY,NP ) 
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300 CONTINUE 
310 IBTA=1 


STORE BETA TEMPORARILY ON DISK SO WE CAN 
COMPUTE THE TRANSPOSE OF TRIX 


REWIND 44 
WRITE! 4 4 ) BETA 
REWIND 44 

STORE TRIX IN BETA 

00 320 M-.l t NP 
DO 320 K=1,NP 
320 BE TA (.M t K ) •= TR I X (.M»'K )' 

NOW TRANSPOSE BETA (.OLD TRIX) 

00 330 M - 1 »• N P 
DO 3 30 K a 1 » NP 
330 TRI-X(.M,K )-BETA(.K,M) 

1 8 T A= 1 

RESTORE BETA 

READ ( 44 ) BETA 

INVERT TRIX 

CALL M I NV (TRIX ,NP ,'AK ,'LL , MM ) 

CALL CAT$W(4,I) 

GO TO ( 340 * 350 ) , 1 

340 WRITE ( IP, 480 ) 

CALL SSS (. I k I X ,NP ) 

350 IF (-NF L A P—1 ) 360,370,370 
360 WRITE ( IP , 490 ) NAME 
WRITE ( IP, 500) 

CALL AAA(.Y,NP) 

WRITE! IP, 510) 

CALL AAA (.ALMAX » NP ). 

SET UP FUR CL MAX LOOK UP 

IY=1 

REYUN=0. 

A L P PI A = 9 9 9 . 

CLL=999. 

KEY N =999. 

LOOK UP CL MAX VALUES 

X M A X = 0 . 

CALL BRIDGE CL MAX, 1,NP» 1, IY, I E , R F Y , TAU , CAMB , AL PH A , 1 , 1 , XMAX ) 
WRITE ( IP, 520) 

CALL AAA (CL MAX , NP ) 

370 IF(A'PLAP) 3 80,380,390 
380 WRITE! IP, 530) 

CALL AAA (-TAU , NP ) 

WRI IE( IP,. 540) 

CALL aAA(KEY,NP ) 

WRITE ( IP, 5 50) 

CALL AAA (.CfNP ) 

WRITE ( IP, 560) 

CALL A A A (. E P S , N P ) 

WRITE ( IP, 570) 

CALL AAA (.CAMB, NP ) 

IF (.; IF LAP— I ). 420, 390, 390 


CHECK FUR DUMP 
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300 CALL DAT SW (.3 * I- ) 

GO TO ( A 00 ,'4 20 ) , I 
400 C ALL DATSW (.3 » I ) 

GO TO (-420 , ! 4 10 ) , I 
4.10 WRITE UP, 580) 

STAR=1STAR 
CALL ZZZ(STAR) 

WRI-TE! IP,' 450) 

CALL AAA (.TONY, -NP ) 

WRITE ( I P ,440 ) 

CALL AAA LGENb »NP ) 

WRITE (IP, 460) 

CALL AAA (.HOPP, 'NP) 

420 RETURN 

430 FORMAT UOX,.L 1HMATRIX BETA) 

440 FORMAT (-10X , 4HGENE 1 
450 FORMAT (-10X,'4HT0NY ) 

460 FORMAT ( 10X, 4HL0PP ) 

470 FORMAT L10X,'5HS I GKA ): 

400 FORMAT ( 1 OX , 11 ft MATRIX TRIX) 

490 FORMAT ( 1 H 1/ 1H0/ 35X , 2 5 A.2 / 1 X ). 

500 FORMAT MX/10 X, 1 17 FI SPA NW I- SE STATIONS) 

510 FORMATUOX, 9HALPHA MAX)! 

520 FORMAT ( 10X,7H CL MAX) 

530 FORMAT ( 10X ,'30H TH I CKNESS / CHORD DISTRIBUTION) 

540 FORM A F( 10X,<24HSECT ION REYNOLDS NUMBERS) 

550 FORMAT ( 1 OX, '18HCH0RU DISTRIBUTION) 

560 FORMAT (.10X, 1 lb H TwIST DISTRIBUTION) 

570 FORMAT ( 10X, 19HCAMBEK DISTRIBUTION) 

580 FORMAT ( 1 0 X , 5 H I S T A R ). 

END 

OVERLAY ( DL I NDA, 2 , 0 ) 

PROGRAM UO 

DIMENSION ARRAY ( 5 , 25 ,• 1 2 ) , C (- 1 9 ) , E P S ( 19 ) ,• TR ANS ( 1 9 ) , RE Y ( 19 ) , fc I A ( 1 9 ) , 

1 HOPP ( 19) * C L M AX (. 1 9 ) ,ZHERL(2,6) ,WHERL(2*6) ,Y(. 19) ,TAU( 19) ,6ETA( 19,19) 
2, TRIX (.19, 19) ,KAZZ(6) ,MZC0L(.6) ,'MAXX(-6) ,MXC0L(6) , C6G ( 19 ) , C VAL ( 19 ) , 
3ALPG( 19 ) »CBC( 19) , A.LPHU (. 1 9 ) , A LPH ( 19) ,ALPHZ( 19) ,ALPHE ( 19) , DELTA* 19) , 
4ALPHV (20 ) , YHERE (2,6 ) , MiYCOLt-fc ) , MAYY (-6 1, CM ( 19 ) , XHERE ( 2 , 6 ) , MA WW ( 6 ) , 
5MWCUL ( 6 ) 

COMMON INNOW , I Si" 1 T (.3 ) , 'ALPHA , REYN ,CLL ,'RE YON ,'XMAX , ALMAX ( 19 > ,CLMAX ,C , 

1 E PS , TRANS, KEY, 'ETA, 'HOPP , Z HERE , WHERE ,'ARP A Y , Y fBET A , TF AC , TR I X , T AU , MAXX 

2 , MAZZ , MXCOL , MZCOL , ASP EC , T APER ,BF ,'K E YND ,'D I SCR , P I E R , CR8 , Q , TSTAX , EDGL 

3 ,'S I G , aL PI'R , NF L AP ,wLVL , NP , I Y , I Z I K , 1P,TG, I ST AR , A , B , H , TAUT , TAUR , 

4T W I ST , R , BWX , YFTHRE ,MYCO.L , MAYY, FLAP, T 0 N Y (. 19), <TWISA,X,Z, CM, ACC, XHERE, 
5MWC0L , M.AwW ,'CAMB ( 19 ) , CAMBR , C AMB T , UUMY 1 , DUM Y2 , NAME ( 25 ) , AHEKE (2,6) , 

6 M A A A ( 6 ) ,M, ACOLt.fi) , EM I ERE (2,6) , MABB (.6 T , MBCOL ( 6 ) , CHERE ( 2 , fa ) ,MACC(6) , 
7MCCUL ( 6 ) , UHL RE (-2 ,6 ) , MADD ( 6 ) , MDCUL ( 6 ) ,'CVAL,ALPHU,CBC,CBG, DELTA, 
8ALPHZ , ALPH, ALPHV, ALPHE, ALP HB , IP RAB ,' I S S ,' I Q T ,'D I R , I AGN 
10 CALL M A I N 2 
CALL MAIN4 
If (- lur-2) 10,20,20 
20 RETURN 
END 

SU BRUu T I NE Mi A 1 N 2 


SUBROUTINE FOR THE CASE WITHOUT A FLAP 

U 1 MENS I UN ARRAY (5, >2 5, 12) »C(.19) , EPS (.19) , TRANS ( 1 9 ) , RE Y ( 1 9 ) ,F T A ( 1 9 ) , 
1H0PP ( 19 ) ,CLMAX (.19 ) , Z HERE ! 2 , 6 ) , WHERE ( 2 , 6 ) , Y (. 1 9 ) , T AU ( 19 ) , BETA ( 19 , 1 9 ) 
2 , TRIX ( 19, 19 ) ,MAZZ (.6 ) ,MZCOL (.6 ) ,'MAXX(.fi ) , MXCOL ( 6 ) , C BG ( 19 ) , C VAL ( 19 ) , 
3ALPG( 19) ,CBC( 19) , ALPHUM9) , ALPH ( 19 ) , ALPHZ ( 19 ) , ALPHE (19),DELTA(19), 
4AI PHV( 20 ) , YHLRFI 2,6 ) , MY COL (.6 ) ,'M A Y Y ( 6 ) ,CM ( 19 ) , XHERE ( 2 ,'6 ) , MA W W ( 6 ) , 
5MWC0L ( 6 ) 

COMMON INfJDrt, ISv.lTt.3) , ALPHA , R F YN ,C L L ,'RF YON rXMAX , ALMAX ( 19) ,CLMAX,C, 
1LMS, TRANS, KEY, ETA ,HUPP, Z HERE , WHL RE , 'ARRAY , Y ,<BE TA , TFAC ,-TR I X , TAU, MAXX 

2 , MAZZ, MXCOL , MZCUL > ASPLC , T APER , BF ,'RE YND ,'D I SCR , P I FR ,CRB , Q , TSTAX , EDGE 

3 , S I G , ALPHR ,'NF LAP ,'NLVL , NP , 1Y,IZ,IR,IP,1G, i ST AR , A , B , H , T AUT , TAUR , 
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A TWISl , R , BW X , YHERE , MYCO.L » MA Y Y , f'L AP , TONY (, 1 9 ) ,XW I SA , X , Z , CM , ACC , XHERE , 
9 HW C U L , M A W W , C A M B ( 19 ) , CAMBR ,\CAMBT , DUM Y 1 , DUM Y2 , NAME ( 25 ) ,-AHERE (2,6) » 

6 h A A A { o ) , MACOL («o ) , BH E K E ( 2 , 6 ). , MAli B ( 6 ) , MBCO L ( 6 ) , CHERE ( 2 , 6 ) *MACC(6) , 
7MCCfJL(6) f LihhRfc (2,6) , MADD ( 6 ) , MDCOL ( 6 ) ,*C VAL f ALPHU ,CBC f CBG , DELTA » 
BAL’PHZ , AL PH ,< U PHV , A.LPHL ,ALPHB , I-PR AB ,-I S S rl OT , D I R , IAGN 

LOOPING PARAMETER. IE IR (THE LOGICAL UNIT 
NUMBER OH THE READER) IS SET PROPERLY AT 
1 1 * I S TIME THEN THIS IS THE FIRST ENTRY INTO 
THIS SUBROUTINE. IT IR IS 100 THEN WE WISH 
TO ITERATE AGAIN TO HIND CONVERGENCE. 

IF UR-8) 20,10,60 
LO I W“ 1 

I P R A B = I 
I S= 1 
ISSUS 

INPUT OF ALPHA VALUES TO BE USED FOR THIS 
RUN. A LIST OF ALPHA VALUES MAY BE 
TERMINATED IN TWO WAYS. (1) INCLUDE ALPHA 
VALUES LARGE ENOUGH TO CAUSE A STALL. (2) 
END LIST OF VALUES WITH A 99. 

READ( I K f 3 2 0 ) ALPHV 
I R - 1 0 0 
GO TO 30 

IF IR IS ZERO THEN WE WISH TO TAKE THE NEXT 
ALPHA VALUE 

20 IS=ISS 

I PR Ab~ I PRAB + i 
30 D I R = G • 

IQ r = i 

ALPHB=ALPHV ( I PRAB ) 

TEST FOR END OF ALPHA VALUES 

IF LALPHB-99. ) 5u, 40, 90 
AO CALL FXIT 
50 kRI-TM IP, -330) NAME 

WRITEt IP, 3 AO) ALPHB,DISCR, A,B, ASPEC,H, ALPHR, TAUT 

WRI rt (.IPf 350) TAURfTWIST, RfTWISAfBF, TAPFRfFLAP, REYND, Xf Z 

WRITE! IP, 360) 

I AGN=0 
60 I T R = 0 
IS=IS S 
I Y= 1 
L OC LR = 2 
DO 70 K = 1 , N P 

70 A L P G ( K ) = ALPHB+.AL PHR + bPS (K) +.ALPHB*TFAC* (-TRANS ( K ) - 1 . ) 

SWITCH NUMBER 3 IS USED FOR AN INTERNAL 
DUMP OF ARRAYS CUMPUTED DURING ITERATION 
PROCESS 

CALL L AT SWi (.3 , JUNK ) 

GO TO ( 80,90) ,• JUNK 
80 WRITF( IP, 370) 

CALL AAA {.AL PG ? NP ). 

FIND C L ( K ). 

90 CLL =999 . 

REY0N=999. 

XMA X = 0 . 

CALL BRIDGLCVAL, I ,'NP , 1 , I Y ,■ 1 E , RL Y ,<T AU , C AMB , A L PG , 1 ,NP ,CLL ) 

IT TIL) 100,100,310 
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100 DC 110 K = 1 f'NP 
A K = K 

110 CBG(K)=CVAL<K).*( ASPQC/ (.ASP EC+ 1 • 6 ) ) *C LK ).*CRB* ( 0 . 5+ ( 1 . +TAPER ) *S I N < AK 
1 * P I- C R ) / ( 3 • 1 A 3 5 9 * C (*K ) ) ) 

CHECK FOR DUMP 

CALL DATS* (-3, JUNK). 

GQ TO M20»-130). f JUNK 
120 WRI-TU IP, 380) 

CALL AAA(.CVALrNP) 

LOOK UP ALPHA VALUE FUR ZERO LIFT 

130 L0CER=3 

ALPHA-999 • 

CLL-.OOG1 

REY0N=999. 

X M A X = 0 • 

CALL bPIDC( ALPHZ.-1 t-.MP»a f 1 Y ,• 1 1 ,< Rfc Y , T AU , C AMB ,'ALPHA , 1 , 1 ,iALPHA ) 

IF ME ) 140, *140,310 
1A0 DC 160 K= 1 1 NP 
S I G=0 , 

DO 190 M=1,NP 

150 S I C=CbG ( M ) *HET A ( M ,«K ) + S IG 

A L P Hu OK ) = S I C * (. 1 . +. I f- AC * {. T R A N S ( K ) - 1 . )! ) 

ALPH(K) =ALPG( K ) -ALPHULK ) 

160 ALPHE<K) = (.ALPH(.K)-ALPHZ(K):*{ 1,-EDGE) ) / E DGE 

LOOK UP CL VALUES FOR A LIST OF ALPHA 
VALUES 

LQ.CER = 4 
CLL=999. 

REYUN=999. 

X M A X = 0 • 

CALL LkICGLCVAL, 1,NP, I, I Y ,- 1-E , KEY ,-TAU ,«C AMB , ALPHE , 1,NP,CLL) 

IF ML ) 170,' 170, 310 

170 DO 180 K - 1 ,»i\P 

CbC ( K ) =C V A L ( K )>C (.K ) *CRB/KWX 
180 D F L T A ( K ) = C B C (K)-CDG(K) 

CHECK FOR DUMP 

CALL UATSWL3, JUNK) 

GO TO { 190,200), JUNK 
190 WR 1* T E ( I P , 390 ) 

CALL A A A (-C B G , N P ) 

WR FTL {. I P *-4 1 0 ) 

CALL AAA (, A L P HU ,«N P 1 
w R I- r E ( 1 P , 4 2 0 ) 

CALL AAA (, ALP HE , IMP ) 

WR ME ( I P, 3M0 ) 

CALL AAA(XVALt-F\}P) 

WRIIEI IP ,420) 

CALL AAAtCBCtNP) 

WRI TE UP,430) 

CALL AAALDELT A,-NP ) 

CHECK TOLERANCE 

200 DO 220 K=1,-NP 

IF(ABS<DELTA(K) )-DISCR) 210,210,230 
210 1F(K-NP) 220,300,300 
220 CONTINUE 
230 DO 250 1=1, NP 
SUM=C. 

DO 240 J - 1 ,■ N P 

240 SUM = TRIXM ,' J ) * D h L T A ( J ) + SUM 
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250 CBG (. I )=CBGU )+.SUM 


CHECK FOR DUMP 



CALL DATSWL3, JUNK) 

GO TO (.260, -270), JUNK 
260 WRI-TE( IP, -390) 

CALL AAA (.CBG , NP ) 


REPEAT CYCLE 


270 ITR-ITR+1 

IF (. I TR— 30 ) 130,130,280 


IF UNABLE TO CONVERGE AFTER 30 ITERATIONS 
DUMP DELTA VALUES, C VALUES, AND TABLE 
PRESENTLY IN CORE bElNG USED FOR LOOK UP 


280 WRITE! IP, 440) 

WRITE! I P *’4 30 ) 

CALL AAA (.CELT A, NP ) 

WR I TE ( I P ,380 ) 

CALL AAA (.OVAL, NP) 

DO 29C J= 1 , NLVL 
N R ■= M A X X ( J) 

NC=MXCOL U ) 

TAUX=wHERE (.2 , J ) 

WRITE! IP, 450 ) J , J ,NK ,-NC , TAUX 
DO 290 K = I ,'NR 

290 WRITE! IP, 460) {.ARRAY (J , K , L ), L= 1 , NC ) 

CALL EXIT 

300 WRI-TE!. IP, 470) I- T R , A L P H 8 
ISS=IS 
RETURN 

3 LO WRITE! IP, 480) lE.LOCLR 
CALL EXIT 

320 FORMAT! 1CF8..0). 

3 30 FORMAT! 1 1! 1/ 1 IiO / 35 X , 23 A2 / IX T 

340 FORMAT! IX, 40! 3H./.):/iH0/15X,3UH BODY. ANGLE OF ATTACK, DEG. .. =,FLO. 

1 2 , 1 0 X , 3 0 H V A L U E OF DISCRIMINANT. . .! .. = , F 1 0 . 6/1 5X , 30HB0D Y HEIGHT / 

2 SPAN = , F 1 0 . 2 , 10X,30HB0DY. WIDTH / SPAN ‘ =,F10. 

32/16X, 30 F. ASPECT RATIO -,F10.2,10X,30HWING HEIGHT / 

4 SPAN =,F10.2/15X,3UHWING BODY INCIDENCE, DEG . . =,F10. 

5 2 *• 1 0 X , 3 0 H i IP THICKNESS CHORD.. .... =,F10.2) 

350 F0P.MA1 ( 15X, •30HK0GI THICKNESS CHORD .... = , F 1 0 . 2 , 10X , 30HGE0ME TR I 

IC TWIST, DEG .... =,'F10.2/15X,'30HNUMBER OF SPANWISE STATIONS. 

2.F10.2, 10X, 30HAER0DYNAMIC TwIST, DEG . . . = , F 10 . 2/ 15X , 30HFLAP SPA 

3N /. „ 1 NG SPAN. . . . =, F10.2, 10X,30FTAPER RATIO = 

4,F10.2/15X,30FFLAP SETTING, DEG = , F 1 0 . 2 , 1CX,30HREY MOLDS 

5 NUMitR =,! ; 10. 2/15X, 37 H COORDINATES OF MOMENT REFERENT 

6 E POINT , 14X,2Fa=,F10.2,9X,2FZ=»F10.2) 

360 FORMA K IHO/ LX , 40 !.3H . / .. ) / IX ) 

370 FUrMAT ( IC'X , 4FALPG ). 

380 FORMAT !.10X,’4HCVAL) 

390 FCRMA1 (.1UX,3HCBG) 

400 FORMAT ( luX , 5HALPHU ) 

410 FORMAT! 10X vSHALPHE ) 

420 FORK m T (■ 1 OX , 3HCBC ) 

430 FORT A I ( 1 0 X ,■ 6 F D E L T A ) 

440 FORMAT (.IX, 48IUJNABLE TO CONVERGE AFTER 30 ITERATIONS ABURTED) 

450 FORMAT! 1 HO/ IHO, 7H ARRAY ( , I 2 , 1 1H 1, IT TO ( ,'I 2 , 1H , , I 2 , 1H , , I 2 , IH ) , F 1 

10. 5/1 HO ) 

460 FORMA I ( 12F 10.4) 

470 FORMAT (. 1 0 X *• I 2 , 6 3 H ITERATIONS REQUIRED TU CONVERGE FOR ANGLE OF AT 

ITACK EQUAL TO , FB.l/IX). 

4F0 FORMAT MX, 10HERR0R CODE , I 2 , 1 X ,■ 10HA T SEC T I ON , I 3 , 1 X , 3 2H I N PROGRAM, i 

1XECIJI ION TERMINATED) 

E N D 
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SUBROUTINE MAIN4 


CONTINUATION OF SUBROUTINE MAIN2 

DIMENSION ARRAY! 5, 25, 12) ,C (.19 )., EPS 1 19 ) ,»TRANSC19) ,REY( 19) , E T A ( 19) , 
1H0PP( 19) ,-CLMAX (.19 ).,-ZHEKE ( 2,6 ) ,< WHERE ( 2 , 6 ) , Y {. 1 9 ) » TAUC 19) , BETA ( 19 ,‘19 > 
2,XRIX( L9 , L9 ) , MAZZC6 ) ,.MZC0L (.6) ,«MAXX(.6) ,MXCOL (6) , ALPG( 19) ,ALPHE ( 19) , 
3 CL M9 ) f-CDU (.19 ):,CDOC( 1 9 ) , ALPHU (. 19 ) ,-CBG( 19 ) , CBG( 19) , DELTA! 19 ) , ALPHZ ( 
A 19) , ALRHU9) ,ALPHV(-20) ,XVAL( 19) , YHER E ! 2 , 6 ) ,*M YC OL ( 6 ) , MAYY ! 6 ) ,CM( 19) 

5 ,'XHERE ( 2 , ! 6 )! , MAWW (.6 ) , MwGOL L6 ) 

COMMON INNOW, IS WIT (.3) , A LPHA , RE YN , CL L ,'R E YON , XMAX , ALMAX ( 19) ,CLMAX,C, 
IE PS ,XRANS,REY,ETA,HOPPfZHERE, WHERE, ARRAY, Y, (BETA, TF AC *TR IX, TAU,MAXX 
2,MAZZ,MXCOL,MZCOL, ASPEC, TAPER „BF ,-REYND ,*D I SCR , P I ER,CRB ,Q, TS TAX, EDGE 
3 ,*S I G , ALPHR ,‘NFLAP t'N.LVL ,'NP , IY, IZyIR, IP,*IG, I ST AR , A , B , H , T AUT , TAUR, 

4 T.W I S T »R,BWX,YHERE , M YCO.L , MA YY , FL AP , TONY U 9 ) ,<TW I SA , X , Z , CM , ACC , XHERE , 
SMWCOLtMAWWfCAMB! 19) , C AMUR ,< C AhbT , DUM Y 1 , DUM Y2 , NAME ( 2 5 ) ,AHERE(2,6) , 

6 M A A A ( 6 ) , MACOL (,6 ) ,-BHERE 02 , 6 ) , MABB {,6 Y , MBCO L ( 6 ) ,CHERE(2,6) ,MACC(6)‘ , 

7 M C C 0 L (6 ) , DHERE (-2 ,'6 ) , M ADD ( 6 )■ , M DC OL ( 6 ) ,‘C VA L , A L PHU , CBC , CBG , DELT A , 

8 ALPHZ , ALPH, ALPHV ,-ALPHE , ALPHB, 1-PRAB,-I S S „ I Q 1* ,*D I R t I AGN 
IS-ISS 

I K T = 2 

DO 10 K = 1 , N P 
10 CL(K)=CVAL(.K) 

I Z = 1 
LOG E R= 5 

COMPUTE PROFILE DRAG COEFFICIENTS 

C L L =999 . 

REY0N=999. 

X M A X = G . 

CALL BRIDG (.CDO , 1 , NP , 2 , I-Z , I fc , RE Y , T AU , C AMB , CL , 1 ,NP ,CLL ) 

DO 20 K=1,NP 
20 CL(.K)=CVAL(.K) 

COMPUTE QUARTER CHORD PITCHING MOMENT 
COEFFICIENTS 

C L L = 999 • 

R E Y C N = 9 9 9 » 

XMA X = 0 • 

CALL BRIDGLCMf 1,NP, 3, IW, lEfREY, TAU,'CAMBtCL,'l,NP,CLL) 

DO 90 K ~ 1 , N P 

CHECK IF SECTION STALLED 

IF ( ALPHc (,K ) - A L M A X ( K ) ) 40,30,30 
30 WRI TL! IP, 260) K , AL P HE !.K ), AL MAX (.K ) 

I S = 2 
I S S = I s 
IKT = 1 

COMPUTE SECTION PITCHING MOMENT 

40 SUM 1 = CUS (-3 . 1 4 LO 9 / 1 80 .* CAL PHB-AL P HU !.K ) ) ) 

SUM2 = 5 I N (.3 • 14 159/lfaO.*!. ALPHB- ALP HULK ).) ). 

CM(K)=CM(.K)-X/C(K)^(CL(.K)*SUMl + CDa(.K)*SUM2)-Z/C(K)*(CL(K)*SUM2-CDU 
1 ( K ) * S U M 1 ) 

5 0 C DOC ( K } =C DO ( K ) *C (.K )>CKB/BWX 
WRI TE < IP, 270) 

CALL A A A (X M , N P ) 

CHtCK FUR DUMP 

CALL UATSW(.3 f JUNK) 

GO TO (.60,70),- JUNK 
60 WP I TE ( I-P , 280 ) 

CALL AAALCLMAXf NP) 

WRITE! IP, 290) 


175 


non noon 


CALL A.AA(XDOC,NP) 

70 WRl-TE( IP f 300 ) 

CALL A A A (, A L P H E , N P } 

WRHb< IP, 3 10) 

CALL AAA tCPn f NP J 
DU 80 K = 1 1 N P 

8 0 AL PG ( K ) = ALP HU (.K : ) *C L ( K ) *3 • 14159/180. 

WRI-Tfc( I P, 320) 

CALL AAA(.ALPG,NP ). 

WRITE ( IP, 330) 

CALL AAAtCLfNP) 

COMPUTE OVERALL LIFT, DRAG, AND PITCHING 
MOMENT COEFFICIENTS 

SUM 1=0. 

SUM 2 = 0 - 
S U Mi 3 = 0 • 

S U ft 4 = 0. 

DO 90 I=I,NP 
SIJM1=C8CU ). * E T A ( IJ+SUM1 
SUM2 = CBCU )> ALP HUM ) *E T A ( I )+SUM2 
SUM3=.C DOC ( I ) * E T A ( I )+SUM3 
90 S U M 4 = C M { I ) *C( I ) **2*ETA(. I )+SuM4 
Q=ASP EC*BWX * *2 
COL I F r = 0*SUMl 
001=0. 017453*0*. SUM 2 
CDP=Q*SUM3 
CD = CD.I + CDP 

Z M = AS P E C * BW X * CRD / ACC * SUM A 

W R I- T E ( IP, -34 0 ) ALPHBfCDI-fCCLIFT ,CDP ,-ZM,CD 
GO 10 (.2 50,-100), IS 

DEFINE EXACT STALL ANGLE OF ATTACK 

100 I AGN = I AGN+ I 

I F ( IAGN-2) 120, 120, 110 
110 IfCOlK) 130,120,140 
120 GO rU ( 160,180) , IKT 
130 GO TO (160,150) , IKT 
140 GO TO ( 150,180) , IKT 
150 1 0 T =2 
160 I Y= 1 

DU 170 K - 1 , N P 

170 I ON Y ( K ) = CLM A X (.K ) - C L ( K ) 

WRITE ( IP, 350) 

CALL AA A (. TON Y , NP ) 

GO TO (180, *240, IQT 

180 GO TO ( 190, 210, 210, '2 10, 240 1, I AGN 
190 I PR AB= I PR AB- 1 

IF (.IPRAb-T); 200,2 10,200 

200 AL PHb= ALPHB- ( ULPHVt I PR AH »A L PHV U PR AB- 1 ) )/2. ) 

AL’PHb=ALPHB + 0.001 
GO TO 220 

2 LO AL PHB= IALPHB+ (.- 1 . )>.* IKT *0.2 ) 

A L PHB = A L PHB + 0*001 
D I R - ( - 1 • ) **IKT 
220 WRITE ( IP, 370) NAME 
WRITEUP, 360) ALPHB 
I PR AB = I PR A B + I 
I R= 100 
230 RETURN 
240 IQT =2 
RETURN 
250 I R-0 

RETURN 

260 FORMAT MOX, 11 H STALLED AT f *I-3 f 32H POINT SECTION ANGLE OF ATTACK =, 
1F10.3, 10X,30HMAX* SECTION ANGLE OF ATTACK * t F10.3) 
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270 FORMAT (.10X ,'35HSEC T I UN PITCHING FOMENT COEFFICIENT) 

280 FT! KM. AT ( 10X t'SHCLMAX ) 

290 FORM. A 1 ( 10X,'4HC0UC ) 

300 FUKMAT(.1X/10X,33HLFFLCTIVE SECT ION ANGLE OF ATTACK) 

310 FORMAT' ( 10X , 32HSEC T I ON PROFILE DRAG COEFFICIENT) 

320 FORMAT! lOX , 32HSECT ION INDUCED DRAG COEFFICIENT) 

330 FORMAT (.1 OX ,'40HD 1 S TK I BUT I UN OF SECTION LIFT COEFFICIENT) 

340 FORMAT (1 X » 1 20 (. LH. ) / 1H0 / 1 5X »'30Hi30DY ANGLE OF ATTACK, DEG. . =,F10.2 
1, 1GX, 30HINDUCED DRAG COEFFICIENT . . = ,-F 10.6/ 1 5X , 3CHL I F T COEFFICIE 

2 NT =,‘F10.6,'10X,30IIPK0FILE DRAG COEFFICIENT . . =,F10.6 

3/15X,30HPF I CHI-NG MOMENT COEFFICIENT. = ,'F 10.6, 10X , 30 H TOTAL DRAG COE 
4FFICIFNT . . . = ,' F 1 0 . 6 /.'I H 0 / 1 X »■ 1 2 0 ( 1 H . ) ). 

350 FORMAT ( lOX , 25HSTALL MARGIN U I STR 1 BU T I ON/ IX ). 

360 FORMAT! IX, 40(3H./. 1/ 1X/.10X , 28HADJUSTED ANGLE OF ATTACK T0,F8.2,32H 
1, SEARCHING FOR EXACT STALL POINT ,/lX/lX»40(.3H.Y. )/lX) 

3 70 FORMAT (. 1 H 1/ 1 H0/ 35X , 2 5 A2/ 1 X I 
END 

OVERLAY ( B.L I NO A ,3,0) 

PRUGKAM THKLL 
CALL MAIN3 
RETURN 
END 

SUBROUTINE MAIN 1 


SUBROUTINE FOR THE CASE WITH FLAPS 

DIMENSION ARRAY! 5, 25, 12) ,C (.19) , EPS (.19) , TRANS! 1 9 ) ,KE Y ( 1 9 ) , ET A ( 19 ) , 
1H0PTM 19) , CLMAX (. 19 ) ,ZHCkfc(2,6) ,'WHERE ( 2 , 6 ) , Y ( 1 9 ) , T AU ( 19 ) , BE T A ( 1 9 , 1 9 ) 
2 , TRIX! L9, 19 ) ,MAZZ (6) ,'MZCOL (.6 ) ,MAXX (.6 1, MXCOL ( 6) , XHERE ( 2 , 6 ) , MWC0L(6) 
3,MAWh(6) ,CM( 19) ,CBG( 19) ,CVAL( 19) ,ALPG( l 9 ) , CBC ( 1 9 ) , ALPHU ( 19 > , ALPH ( l 
49) , ALPFZ (.19 ) , ALPHE (.19) , DELTA! 19 ) ,'CL ( 19)., ALPHV! 20) , YHERE ( 2,6) ,MYC0L 
5(6) , MAY Y (.6 ) , CLAUD M9) , CLDE L ( 19 ) , C L AD2 ( 1 9 ) ,C L AD 1 ( 19 ) , F ( 19 ) , ALPC ( 19 ) 
6, TUNY (19) 

COM PUN I NNOW , i Sw I 1 (.3 ) .ALPHA, REYN , C L L , R E YCN , XMA X , ALMAX ( 1 9 ) , C LMAX , C , 
1EPS, TRANS, RE Y, ETA, HUPP ,ZhEKE, WHERE, ARRAY, Y.'BFTA, TFAC,'TR IX, TAU.MAXX 

2 » / '/\ Zz. , I X CO L , NizCOL , ASPEC , TAPER ,'BF , RFYNO ,'D I SCR , P I ER , CRB , Q , TS T AX , EDGE 

3 , S FG , AL PFiR , kF LAP , NL VL , NP, I Y. , I Z , I R , I P ,' I S I S , I S T AR , A , B , H , TAUT , TAUR , 

4 IWI-SI , R, BWX, YHERE, MYCOL,MAYY, FLAP, TUN Y,TW I SA,X,Z,CM, ACC, XHERE, 

5 M w C 0 L , M A W W * C A M B ( IV) , CAMBK , CAMB T , DUMY 1 , DUMY2 , NAME ( 25 > , AHERE (2,6) , 

6 M A A A ( 6 ) , M A C 0 L (. 6 ) , B HERE (. 2,6) , M A B B ( 6 ) ,MBC0L(6) ,CHERE(2,6) ,MACC(6) , 
7MCCOU6) , OH ERE (2 ,'6 ) , MAOO ( 6 )., MDCOL ( 6 ) ,'ALPHE ,'ALPHV , CLAUD , CLDEL , CLAD2 
8 , C.L AD 1 , F , IK I , F'F , LUCtR ,'I PRAB , IK T , IUT , I AGN , I Tl< , I S 

ALPHC=n. 

CL ( 1 ) =0. 

J P-R/2 . 

1 PRAB=U 
l« 1 = K 
I S = 1 


INPUT ALPHA VALUES 


KEAU< IK, 700 ) ALPHV 
I PR AB= I P R A 8 + 1 
D I R=G . 

I QT= 1 
I AGN=0 
I TR = 0 
I Y= 1 


PUT CUBE 1 IN CORE 


CALL UAGETT ARRAY »'l ,'IY ) 
L 0 C E R = 2 
DO 10 K= 1 , N P 
10 AL PG ( K ) = 3 . 

CALL DATSw (.3, JUNK > 

GO TO (20,30), JUNK 
20 WRI TEUP, 710) 

CALL A A A (. A L P G » N P ) 
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LOOK UP CL VALUES 


30 DO 60 K=1,NP 
CLL=999 . 

ALPHA=ALPG(.K) 

REYN=REY(.K) 

T AUX = T AU ( K ). 

REY0N=999. 

XPAX=0. 

CALL ARC (.ARRAY, TAUX , MA.XX , MXCOL , I E , WHERE , NLVL ) 

CHECK FOR ERROR STOP 

IF(IE) 50, 50, AO 
AO WRI-rE< IP, 720) IL.LOCER 
CALL EXIT 
•50 CVAL ( K ) =CLL 
AK = K 

60 CBG(K)=CVAL(K)*(ASPEC/(.ASPEC + l.fa) )*C(.K)*CRB*(0.5+( 1 . + TAPER ) * S I N ( AK 
1*PICR)V (3.1A159*C(K) ) ). 

CALL CATSW1.3, JUNK). 

GO TU (70,80), JUNK 
7.0 WR I T E ( IP, '7 30) 

CALL AAA (.CV AL , NP ) 

80 DO 120 K = L ,• N P 
S I G=0 . 

DO 90 P = 1 , N P 

90 S I G=CBG ( h ) «.BETA ( M ,-K ) +SI G 

ALPHU ( K ) = S 1 G* (.L.+.T FAC*. (.TRANS (K )-l. )'. ) 

ALPH(K) = ALPG(K) -ALPHU (.K) 

LOOK UP ALPHA FOR ZERO LIFT 

LOC ER= 3 
ALPHA=999. 

C LL = 0 . 000 1 
REYN=REY(.K) 

T A U X = T A !J (• K )' 

REY0IM = 999. 

X M A X = 0 . 

CALL ARC (.ARRAY, TAUX, KAXX, MXCOL , I E » WHERE * NL V L ) 

IF ME) 100, 100, AO 
100 ALPFiZ ( K) = A L P F A 

ALPHE(K) = (.ALPH(K)-ALPHZ(K)*( 1 . -EDGE ) )/EDGE 

LOOK UP CL VALUES 

LOCER=A 
C L L = 999 . 

ALPhA = ALPHE ( K ). 

RE YN = RLY (.K ) 

T A U X = T A U (K ) 

REY0.N=999 . 

X P A X = 0 . 

CALL ARC (.ARRAY, TAUX, WAX X,MXCUL , I E , WHERE , NLVL ) 

IF LIE) 110, 110, AO 
110 CVAL ( K ) = CLI. 

C«C ( K ) =CVAL ( K ) ( K ).*CKB/BwX 

120 UE.L TA ( K ) =CBC ( K ) -CBG ( K ) 

CALL DATSw (.3 , JUNK ) 

GO TO (130,'1A0), JUNK 
130 wRI-TEt IP, 7 AO) 

CALL AAA(.CBG,NP) 

W R I T t ( IP, 750 ) 

CALL AAA(.ALPHU,.NP > 

WR ME ( I P , 760 ) 

CALL A A A (. A L P H E , N P ) 

WR Ht( I P , 730 ) 


178 



non non nnn non 


CALL A A A (C V A L , N P ) 

WR HE ( IP, 770) 

CALL A A A (. C B C * N P ) 

WRI-TH I- P 7 5 0 ) 

CALL AAA (.DEL I A ,'NP ) 

1 AO DO 160 K-l ,'NP 

IF(.ABS(DE.LTA(K) )-DI-SCR). 150,150, 170 
150 IF(.K-NP) 160,240, 240 
160 CONTINUE 
170 DO 190 1 = 1, 'NP 
SUM=0. 

DO 160 J = 1 ,'NP 

180 SUM= I R I-X (. I- »' J ) *.DELTA( J ) +.SUM 
190 CBG ( I ) = CI3G (. I J+SUM 
CALL UATSWI.3, JUNK) 

GO TO (.200, '210), JUNK 
200 WRITE ( IP, '740 > 

CAUL AAA (.COG, NP) 

210 l T R= I TR+ 1 

ITERATION CONTROL 

I F (■ l TR— 3 0 ) 80,'60,'220 

IF UNABLE TO CONVERGE, DUMP- ALL 

220 WRITE ( IP, 790) 

WR 1TE ( I-P ,'780 ) 

CALL AAA (.DELTA, NP) 

WR I-TE ( I P , 730 ) 

CALL AAA (.C VAL , NP ) 

DO 230 J-L,'NLVL 
NR=MAXX ( J ) 

NC = MXCOL (.J ) 
r AUX = WHERE (.2 , J ) 

WRITE ( IP, 8 00) J , J ,-NR ,'NC , T AUX 
DO 230 K-l ,'NR 

230 WR I-TE ( I P , 6 1 0 ) (AKRAY(.J,K,L) ,L=1,NC) 

CALL EXIT 

CONVERGENCE 

240 WRTIE( IP, 820) NAME 

WR I-TE HP, 830) I TR ,'ALPHB 

I K T = 2 

SUF=U. 

DC 250 I = L ,'NP 
250 SUM = CiiC ( I ) *,E T A (. I ) +SUM 
Q=ASPEC*6wX**2 

LIFT COEFFICIENT 

COLIFT=Q*SUM 
DO 266 K-l ,-NP 

260 CLADU(K)=CVALl.K) /COLIFT 
CLL=0. 

ALPHA=999 . 

REYN=RfcY(. ISTAR) 

T AO X= T AU (.ISTAR ) 

R F Y C N = 9 9 9 . 

CA.LL ARC (ARRAY, T AUX , MAX X ,MXCUL , I E, WHERE ,NLVL ) 

A 1 = ALPHA 
AL PHA = 999 . 

CALL DAGET (.ARRAY ,'4, 1W ) 

CA.LL ARC (.ARRAY, T AUX , MAY Y , MYCOL , I E, YHERE , NLVL ) 

A 2 = A L P 1 1 A 
SL-E L T = A 1-A? 

ALPG{ JP ) =0 . 

ALPUA = ALPG(.JP ) 


179 



CLL-999. 

RE YN=KEY (-JP ) 

T A U X = T A U (• J P ) 

C 

C LOOK UP f-LAP LIFT (4) 

C 

REY0N=999. 

CALL DAGfcT (.ARRAY , 4 ,1 w ) 

CALL ARC (.ARRAY , I AUX, MAWW , MWCOL , I E , -XHERE , NLVL ) 

CVAL ( JP ) =CLL 
DU 290 K-1,JP 
AK = K 

S TAR= I ST AR 

IF (.K- IS TAR I! 280 » 280 * 2 70 

2 70 C8G(K ) =0.5*CVAL ( JP )>C(K ) *CRB/BWX»( 1 . + SQR1 ( l.-(COS( AK*PIER) /COS( 
l S T A R * P I r R ) ).**2) ) 

GO TO 290 

28 0 CBG( K ) = 0. 5*CVAL ( JP ) *C( K ) *CRB/I3WX* ( 1 .-SQR I ( 1 .-( ( 1 .-COS ( AK*PI EK ) ) / 
l( l.-CUS(SFAR*PIER) ). )**2) ) 

290 CONTI NUE 

DO 800 K - 1 ,• J P 
IR (■ = R 
I R I-= I R I--K 

300 COG ( I R I ) - C W G ( K ) 

IR I=R 

DO 3 I C K = l » N P 
310 A L P G ( K ) = 0 • 

CALL DA F Sin (.3 , I ) 

GO TO (320,330), I 
320 WRI-Ttl IP, 840) 

CALL AAA (.CL ADD, NP ) 

Vv K I T t ( I- P , 8 5 0 ) 

CALL Z Z Z (.CL L ) 
wR I IL( I P , b6 0 ) 

CALL AAA(.CBG,NP ) 

330 I T R = 0 
340 DO 400 K = L N P 
SIG=0. 

DO 350 M = L , N P 

350 S I G=C8G ( M ) * BE T A ( M , K ) + S I G 

ALP HU (K)=SIG*(. l. + TFACM. TRANS(K) — 1. ) ) 

ALPH(K)=ALPG(K)-ALPHU (.K) 

ALPC(K) = SDELT»KOPP(.K ) 

alph(k)=alph(k)-alpc(K): 

CLL =0 . 

AL PHA=999 . 

RE YN=REY ( K) 

REY0N=999. 

TAUX = TAU(.K) 

IF (K- I- STAR ): 380,380,360 
360 IF ( K-IRI+.ISTAR) 370,380,380 
370 CALL DAGET(.ARKAY,4,IW) 

CALL ARCURRAY, f AUX , MAWW , MwCOL , I E , 'XHERE , NLVL ) 

ALPHZ(K)=ALPI!A 

ALPHA = ( ALPH ( K ).— ALPHZ ( K ) * ( 1 . -EDGE ) ) /EDGE 
CLL=999 . 

CALL DAGE I (.ARRAY , 4 ,■ I W ). 

CALL ARC (.ARRAY, TAUX, MAWW, MWCOL, It, XHERt, NLVL ) 

CVAL ( K ) =C.LL 
GO TO 390 

380 CALL DAG ET (.ARRAY ,'l , I Y ) 

CALL ARC (.ARRAY, T AUX, MAXX,MXCOL, I E , WHERE , NLVL ) 

ALPHZ(K)-ALPHA 

ALPHA=(.ALPH(K)-ALPHZ(K).*( 1 .-EDGE ) ) /EDGE 
CLL=999. 

CALL OAGET (.ARRAY ,1, IY) 

CALL ARC (.ARRAY., T AUX, MAXX,MXCOL, IE, WHERE, NLVL ) 

CVAL ( K ) =CLL 

390 CBC(K)=CVAL(K)*C(,K)*CRB/BwX 
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400 DELTA ( K > =CBC ( K>. — CLiGt K ). 

CALL DATSW (.3 , JUNK ) 

GO T(J (4 10, >420), JUNK 
4 LO WRITE< IP, 870) 

CALL AAA (.ALPHU » NP ) 

WRITE! IP,' 800) 

CALL A A A (. A L P C , N P ) 

WRI F E ( IP, 890) 

CALL A A A (• A L P H , N P ) 

write! ip,' 9oo) 

CALL AAAI.ALPHZ ,NP ) 

WRITE ( IP, 9 10) 

CALL AAA (XVAL *NP ) 

WRITE( IP, 920) 

CALL AAA (.CBC ,NP ) 

WRI-TE! IP, 930) 

CALL AAA (.DELTA, NP ) 

420 DO 440 K = 1 » N P 

I F (-ABS (DEL T A ( K ) ) -01 SCR) 4 30,4 30,4 50 
430 IF(K-NP) 440, 500,500 
440 CONTINUE 
450 DO 4 70 I = 1 » N P 
S U M = 0 . 

DO 460 J = l ,'NP 

460 SUM= T K I X (. I- ,• J ) * D E L T A ( J 1+SUK 
470 CBG ( I ) = CBG(.I ) +.SUM 
CALL DArSw(.3, JUNK) 

GO I'Ll (480, '490), JUMK 
460 WRITE (IP, 940) 

CALL AAA (.CBG , NP ) 

490 I T R = I T P. + I 

IF(.irk-JO) 340,540, '220 
500 DU 5 L 0 K- 1 , NP 
510 CLDEL(K)=CVAL(.i\) 

F l = C L 1 1 F L I.JP ) /CLADOUP ) 

F2 = CLUEL(.L)/CLADD(.L) 

DO 520 l< - 3 ,' 1 S T 1 A R 
520 CLAD2(K)-F2*CLADD(.K) 

DO 530 K-ISTAK,JP 
530 CLALU(K)=E1«CLADD(K) 

bl)£LT = CLAUi ( ISTAR) -CL AD2 ( ISTAR) 

DO 540 K=*L,'ISTAR 

540 F ( K ) = ( CLDEL ( K T— CLAD2 ( K ) ) / ll D E L T 

F F= ( CLDE L ( I STAR)— CLADK. ISTAR ) ) / B 0 E L T 
I S T= I STAR + 1 
DO 550 K = I- S T , J P 

55 0 F (.K) = (CLDEL (K)-CLADl(K) ) /ODELT 
DO 560 K = 1 ,■ J P 
1 R I-=K 
I R I- = I R I - K 
560 FUR I > = F- (-K ) 

CALL UATSW< 3, JUNK) 

GO TO (.570, '580 ) , JUNK 
570 WR I TE ( IP, 950) 

CALL AAA (.C LADD, NP) 

WRI-TE! IP, 960) 

CALL AAA (.CLDEL , NP ) 

WRI-TE ( I P ,970 ) 

CALL ZZZ (FI) 

WRI-TE ( IP, 980) 

CALL ZZZ (F 2 ) 

WRITE! IP ,990 ) 

CALL A A A (..C L A D 2 , N P ) 

W R I- T E ( IP, '1000): 

CALL AAA (. C L A D 1 , N P ) 

WRITE! IP, 1010): 

CALL AAA (.F , NP). 

WRITE! IP, 1020) 

CALL ZZZ(FF) 
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580 


590 

600 

610 

620 

630 

64 0 

650 

660 


LOOK UP CL MAX CUfih 1 


CALL DAGET (.ARRAY t 1 »• I Y ). 

RE YN=999 • 

XM AX=0* 

C L L = 999 . 

ALPH A=999 • 

REYON-Rb'Y MSTAR) 

T AUX = TAU (. 1ST A.R ) 

CALL ARC (.ARRAY, T AUX, MAXX ,MXC0L, I E , W HERE , NL VL ) 
CLMNF=XM AX 
XMA X = 0 • 


LOOK UP CL MAX CUBE 4 


CALL OAGET (.ARRAY ,'4 , I W ) 

CALL ARC (.ARRAY, T AUX , MAW W , MWCUL , I E , XH ER E , ML V L ) 
C.LMF = XMAX 
DC.LMA=CLMF— CLMNF 
C LL = 9 9 9 * 

A L P HA= 9 99 • 

REY M = 9 99 . 

X MAX = 0 • 

DO 630 K=l, IMP 
1 R T = R 

IF (.K-ISTAR )! 600,600,590 
I F ( K- I K I + I S T AR ) 610,6 00,600 
X M A X = 0 • 

RE Y0N = RF Y ( K ) 

TAUX=TAU(.K) 


LOOK UP ALPHA MAX CUBE 

CALL DAGET 1 (-ARRAY ,1 ,- 1 Y ) 

CALL ARC (.ARRAY, T AUX , MAX X , MXCUL , I E , WHERE , NLVL > 
GO TO 620 
XMAX = 0 • 

REYON=KFY(K) 

TAUX= f AO (.K ) 


LOOK UP ALPHA MAX CUBE 


CALL DAGET (.ARRAY ,4 , I W ) 

CALL ARC (.ARRAY, T AUX, MAW W,KWCOL, I F , XHFRE , NL V L ) 

F (. K ) = 1 • + F ( K ) *DC LMA/XMAX 

CLMAX ( K } - f (> K ) *. X K A X 

F F = I . + F F * DC L K A / C L M f 

DO 670 K - 1 ,‘N P 

IF (.K—ISTAR). 660,650,640 

I F ( K- I K I + IS T AR ) 660,650,650 

LOOK UP ALPHA MAX CUBE- 


CALL DAGET (.ARRAY,1 , IY) 

ALPH A =999. 

REYN =99 9 • 

CLL=999. 

REY0»M=Rt Y ( K ) 

TAUX-TAU IK ) 

X M A X = 1 0 0 . 

CALL ARC (.ARRAY , T A U X , M A X X , M X C G L , I E, WHERE, NLVL ) 
GO (0 670 


LOOK UP ALPHA MAX CUBE 


CALL DAGET 1 (.ARRAY, 4,1 to) 
XMAX= IOC. 

K E YON = RE Y ( K ) 


(FLAP) 


1 


4 (FLAP) 


1 


4 (FLAP) 
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TAUX = TAU(.K ). 

KEY N =999* 

AL'PHA = 999* 

C LL = 999 • 

CALL ARC LARK A Y» ! T AUX f MAWW » MWCOL » I E t’XHERE , NLVL ) 

670 ALMAX { K ) ^ALPH.Z {.K ) +. {, X M A.X - A L P FI Z U ) ) * EDGE *,F ( K > 

WR IT E ( I P *'820 ) NAME 
•WR I*TE ! IP, 1030); 

CALL AA A (. Y »'NP ). 

WRITE! IP, 1040). 

CALL AAA (.ALMAX , NP ) 

WRI-TE ( I P ,1050 1 
CALL A A A (.CLMAX t NP ). 

WRITE! IP, -LOGO I 
CALL AAA M*AU » NP ) 

WRITE! IP, '1070); 

CALL AAA (.KEY » NP ) 

WRI-FE! IP, L 08 0 ) 

CALL AAA (C ,'NP ), 

WRITE! IP, 1090) 

CALL AAA (-EPS, NP) 

CALL DATSW (.3, JUNK) 

GO 10 (-680, '690), JUNK 
680 WRITE! IP,L100) 

CALL Z Z Z (.DC L M A ) 

WRITE! IP ,10 10); 

CALL AAA(.F,‘NP). 

WRI- 1 E ( IP, 1020 ) 

CALL ZZZ(.FF) 

690 CALL MAIN5 
RE TURN 

700 FORMAT { 10F8..C ). 

710 FORMAT ( IOX, AHALPG) 

720 FORMAT ( IX , 10HERR0R CODE , I 2 , I X 1 OHA T S EC T I ON , I 3 , 1 X , 3 2H I N PROGRAM, t 
1XECUTICN TERMINATED) 

730 FORMAT ( 1 OX , 4HCV AL ) 

74 0 FORMAT (• lOX ,'3HCBG ) 

750 FORMAT ( 10X, 5HALPHU) 

760 FORMAT (• 10X, ; 5HALPHF ). 

770 FORMAT U0X t '3HCBC ) 

780 FORMAT (• 1 0 X , 5 H D (: L T A ) 

790 FORMAT < ? X , 4UHUNABLE TO CONVERGE AFTER 30 ITERATIONS ABORT F D ) 

600 FORMAT MHO/ 1H0, 7H ARRAY ( t I 1 , 1 1H I, 1). TO ( 1 1 * 1H , , I 2 * 1H, t I 2 , 1H ) , F 
110-5/ i H 0 ). 

810 FORMAT (■ I 2F 10. A) 

82 0 FORMAT ( I H 1 / 1 H 0 / 3 5 X , 2 5 A 2 / I X ) 

830 FORMA I ( IX, 33H ITERATIONS REQUIRED TO CONVERGE , 12, 19H FOR ALPUB L'Q 
I V A L 10, F 10. 1/ LHO ) 

840 FORMAT ( 10X, 5FCLADU) 

850 FORMAT M OX f '2HCL) 

860 FOkMAT ! 10X,'3HCBG) 

870 FOKKATUOXt 5HALPHU) 

880 FORMAT (10X , 4HALPC ) 

890 FORMAT {. 1 0 X ,-4 h A L P h ) 

900 FORMAT ( 10X,'5IIALPMZ ) 

910 FORMAT (■ 1 OX , 40C V AL ) 

920 F 0 K M A J ( 1 0 X ,< 3 H C f> C ) 

930 FORMA i !.10X , 5HDFLTA ) 

940 FORMAT! 1 0 X ,• 3 H G B G ) 

950 FORMAT ( L0X,’5HCLADD) 

960 FORMA! ( 1 OX , 5 1 1C L C E L ) 

970 FORMAT (. 1 0 X 2 FI F 1 ) 

980 FORMAT ( 10X,2HF2 ) 

990 F GKF A I (.IOX ,-5HCL AD2 ). 

1000 FORMAT! 1 OX , 5 hCL AO 1 ) 

1010 FOKF AT ( 10X ,-lHF ) 

1020 FORMAT ( 1 0 X , 2 1 iF F ) 

1030 FORMAT ( 1 X / 1 0 X , I 7 M S P A N W I* S L STATIUNS); 
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1040 FORMAT U0X r 10H ALPHA MAX ) ^ 

1050 FORMAT {-1GX »'7H CL MAX) 

1060 FORMAT (10X, 31H THICKNESS / CHORD DISTRIBUTION) 

1070 FORMAT ( 10X..24H SECTION REYNOLDS NUMBER ). 

1080 FORMAT! 10X,'19H CHORD DISTRIBUTION) 

1090 FORMAT ( 1 OX t 1 6fi GEOMETRIC TWIST) 

1100 FORMAT M0X,-6HDUCLMA) 

END 

SUBROUTINE MAIN3 

CONTINUATION OF SUBROUTINE MAIN3 

DIMENSION ARRAY! 5,25 y 12 ) , C (. 19 ) , E PS (, 19 ) , 'TRANS ( 19),REY(19),ETA(19), 
1HOPP ( 19) , C L M 1 A X (. 1 9 ). , Z HERE ( 2. ,'6 ) , 'WHERE (2,6) , Y (.1 9 ) , TAU ( 19 ) , BE TA ( 19 , 19 ) 
2, TR IX( 19,19 ) ,MAZZ (.6 > ,MZCGL(.6 ) , M A XX (.6 T , MXCGL ( 6 ) , CBG ( 19 ) , C V A L ( 19 ) , 
3ALPG! 19) ,C8C( 19) ,ALPHU(.19) , ALPH! 19) ,ALPHZ( 19) , ALPHE! 19) .DELTA! 19) , 
ACL (19) , AL'PHV ( 20 ) ,'XFEKfc (.2 , 6 )., MWCOL (6 ) ,‘MAWW(6) ,CM( 19) ,YHERE(2,6) , 

5 MY C 0 L ( 6 ) ,MAYY(.6) , CLAUD M9) ,CLOEL( 19) ,'CLAD2(19) ,CLA01 ( 19) ,F (19) , 
6ALPC! 19) , TONY ( 19) 

COMMON INNOW, I S a I T (.3 ) ,AL PHA , RE YN , CLL ,'RE YON ,XMAX , ALM AX ( 19) ,CLMAX,C, 
1 EPS, TRANS, RE Y, ET A, HOPP, ZHERE, WHERE, 'ARRAY, Y, 'BETA, TFAC.TR I X,TAU,MAXX 
2»MAZZ,MxCUL , MZCUL , AS PEC, 1 APER , BE , RfcYND ,0 I SCR , P I ER ,CRB ,Q , T STAX, EDGE 
3 , S I C , a L P H R , N F L A P , N L V L , N P , I Y , I Z , I R , I P , I S I S , I S T AR , A , B , H , TAUT , T AUR , 
4TWIST,R,BwX,YHFRF. ,MYCO.L, MAY Y, FLAP, TONY, 'TW ISA, X,Z, CM, ACC, XHERE, 
5MWCUL , FAwW.CAMo! 19) , C AMbR ,-CAM.ti T , DUMiY 1 , UUMY2 , NAME < 25 > , AHERE ( 2,6) , 
6MAAA ( 6 ) , MACUL (.6 ) , BHERE (.2 , 6 ) , MABB (6 ) , MBCOL ! 6 ) , CHE RE ( 2 , 6 ) ,MACC(6) , 
7MCCUL ( 6 ) , DHE Rt (-2 ,6 )' , MADO ( 6 ) , MDCOL ( 6 ) , AL PHE , ALPH V , C LADD , CLDFL , C LAD2 
8.CLAD1 , F,IR I , FF, LOCER, I PR AB, I KT, I OT, IAGN , I TR , IS 
.10 ALPHB=ALPHV ( I PRAB ) 

I F (.ALPHV (. I PR AB ). . EL) . 99 . ) CALL EXIT 

I PRAB= I PRAB+ 1 

WRI TE UP, 560) NAME 

WRITE! I- P ,' 5 7 0 ) ALPHB.DISCR, A,B, ASPEC.H, ALPHR.TAUr 
WRITE! IP, 500) f AUK, TWIST, K, TW ISA, BF, TAPER, FLAP, RE YNO,X,Z 
WRITE! IP, 590) 

20 I K T = 2 

DO 60 K=1,NP 
AK^K 

ALPG ( K ) = ALPHB + A L PHR + EPS ( K ) +.AL PHb*TFAC* !.TR AN S ( K ) - 1 . ) 

IF(.K-ISTAR) 40, 40, '30 
30 I F (. K — I- K I + I- S T A R ) 50, 40, '40 

LOOK UP C V AL CUBE 1 

40 ALPHA = AL PG t.K ) 

CLL=999. 

REYN = REY (.K ) 

TAUX=TAU(.K). 

<EYU:M = 999. 

CALL DAGET (.ARRAY, 1 , IY) 

CALL ARC (.ARRAY, TAUX, MAX X,MXCOL, IE, WHERE, NLVL ) 

CVAL ( K ) =CLL 
GO 10 60 

LOOK UP CVAL CUBE 4 (FLAP) 

50 ALPHA = ALPG(.K ) 

CLL=999. 

RtYN = RLY(.K ) 

TAUX = T AU (.K ) 

REYUN=999. 

CALL DAGLT (.ARRAY ,-4,-IW) 

CALL ARC (.ARR AY, TAUX.MAWW, MWCOL, I E ,XHERE , NLVL ) 

CVAL ( K ) =CLL 

60 CBC (K) = CVAL ( K ) e ( ASPEC / (.ASPEC + 1 .8 ) ) «C(.K )*CRB* ( 0 . 5+ ( 1.+TAPER)*SIN(AK 
1*PI-ER)7(3.14139*C(.K) ) ) 


CHCCK FOR DUMP 
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non nnp non poo 


CALL iMTSw< 3, JUNK) 

GO (U (.70,00),- JUNK 
70 WRI-1 L ( I P , 6 1 0 ) 

CALL AAAMRANS.NP) 

WRITE ( IP, 600) 

CALL ZZZ(.IFAC) 

WR I T E ( I P , 620 ) 

CALL AAA (.AL PG , NP ) 

WRITE ( I P , 630 ) 

CALL AAA (.CBG »f\!P ) 

80 CLSTA=C8G( ISTAR)*hwX/(C( I STAK )>CRB ) 

CHECK FOR DUNP 

CALL UATSV-M.3, JUNK) 

GO 10 (90,100)!, JUNK 
00 WRITE (IP, 640) 

CALL ZZZ(.CLSIA) 

LOOK UP ALPHA CUBE 4 (FLAP) 

100 CALL OAGEI (.ARRAY, 4, IW) 

ALPHA* 999 . 

C L L = 0 . 

REYN = KEY (. I S T AR ) 

TAUX=TAU(. 1ST AR) 

REYUN=999 . 

CALL ARC (.ARRAY, TAUX, toAWW, MWCOL , 1 1 ,-XHtRE , NLVL ) 

A 3 = A L P l-i A 

CALL 0ATSWL3, JUNK) 

GO TO ( 1 10 v 120), JUNK 
110 WRITE! IP, 6i>0) 

CALL 7.171 A3) 

LOOK UP ALPHA CUBE 1 

120 CALL UAGfcT {.ARRAY ,1, IY ) 

CLL=U . 

ALPHA=999. 

CALL ARC (.ARRAY , T AUX , IX AX X , P XCOL , I E , WHERE , NLVL ) 

A4* ALPHA 
CLSTU=CLSTA/FF 
CALL UATSW(.3, JUNK) 

GO TO (130,140), JUNK 
130 WRI T t ( IP, 660) 

CALL ZZZLCLSTU) 

WRITE! IP, 670) 

CALL Z Z Z (• A 4 ) 

140 CALL U A GET (.ARRAY, 4, lw) 

AL PhA = 999. 

CLL=CLS TU 

RE'Yn=RLY USTAk ) 

T AUX = T AU (.1 STAR ). 

RE'Y0N = 999 . 

CALL ARC (-ARRAY. , T AUX , MAWW,MWCUL , I E ,'XHERE , NLVL ) 

ALP HX= ALPHA 

A1=EDGE*FF* ( ALPHX-A3 )+A3 
CALL DATSWL3, JUNK) 

GO 10 (.150, ’160), JUNK 
150 WR 1-Tfc ( IP, 680) 

CALL ZZZ(.Al) 

160 CLS TU = CLST A/F (.I-STAR ) 


LOOK UP ALPHA CUBE 1 


CALL DAGET (.ARRAY, 1 ,1 Y) 
ALPHA=999. 

CLL=CLSTU 
REYN = KFY(. ISTAR) 
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TAUX=7AU(ISTAR) 

REY0N=999 . 

CALL ARC(.ARRAY.,TAUX,MAXX,MXCOL, I L , 'WHERE , NLVL ) 
ALPHX=ALPHA 

A2=EDGE«F ( I STAR ) * (,AL PHX- A4 ) + A4 
CALL DAT SW (-3 » JUNK ) 

GO T 0 ( 170, 180) , JUNK. 

L70 WRI-TE( IP, 690) 

CALL ZZZ(.A2) 

180 SDELT=A2-A1 
I TR = 0 

DO 240 K- 1 ,'NP 
S I G - 0 . 

DO 190 to- L ,'NP 

190 SIG = C8G(M)*BETA(M,K)+SI-G 

ALPHU(.K) =SIG* (.1.+ ri AC« (.TRANS (K)-l. )'. } 

ALPIHK ) =ALPG(K)-ALPHU(.K ) 

ALPC ( K ) = SDE LT*,HOPP (.K ) 

ALPH ( K ) =A.LPH ( K ) - A L P C ( K ) 

C L L = 0 . 

AL PHA=999 . 

KEYN = I<EY (.K ) 

TAUX= r AI.J (-K )! 

REY0,M = 999. 

IF (K-ISTAR): 220,220,200 
200 1F(.K^IKI + ISTAK) 210,220,220 
210 CALL OAGET (.ARRAY , 4 , I Z ) 

CALL ARC (.ARRAY, T AUX, MAW W,MwCOL, I E , XHERE , NLVL ) 

ALPHZ (K ) = ALPHA 

ALP N = ( ALPHIK )-ALPHZ(K )*< l.-F (K ) "EDGE ) ) / ( EDGE*F ( K ) ) 
C L L=999 . 

A L PHA = A L PN 

CALL ARC (.ARRAY, T AUX , MAWW , MwCUL , I E , XHERE ♦ NLVL ) 
CVAL(K)=CLL*F(.K) 

GO 10 230 

220 CALL UAGLTl ARRAY, It IY) 

CALL ARC (.ARRAY, T AUX , MAXX , VXCCL , I E , WHERE , NLVL ) 
ALPHZ(K)=ALPHA 

ALP N = ( ALPH (-K ) - ALPHZ ( K ) * ( l.-F (K ) "EDGE ). ) / ( CDGE*F ( K ) ) 
ALPIiA = ALPN 
C L L = 9 9 9 . 

CALL DAGET (-ARRAY , 1, IY) 

CALL ARC (.ARRAY, I AUX , MAXX , MXCOL , I E , WHERE , NLVL ) 

CVAL ( K ) = C.LL"F (.K ) 

230 CBC < K ) = C V A L ( K ) *C (-K ! "CRB/BWX 
240 DEL TA ( K ) =CBC ( K ) -Cl'G ( K ) 

CALL DATSW(-3, JUNK) 

GO TU (.290 ,260), JU-IK 
2 60 wRI-rt ( IP, Z 2 0 ) 

CALL A A A (. A L P ti U , N P ) 

WRIIL ( IP, ZOO) 

CALL A A A (. A L P 1 1 , N P ) 

WR I I E ( I P , Z 1 0 ) 

CALL AAA (.ALPC,. MP) 

WRME ( IP, 730) 

CALL AAA (.CVAL ,'NP )' 

WRITE (-IP, ZAO ) 

CALL AAA (.DEL I A, -g P V 
260 DO 280 K = L , NP 

IF(ABS(DELTA(K) )-DISCR) 270,270,290 
2 70 If(.K-NP) 280,360,360 
280 CONTINUl 
290 DU 310 I = L , N P 
SUK=0. 

UCJ 300 J = 1 »• N P 

300 SUM=TKIX (1 ,'J ) "DEL f A ( J J-+SUM 
310 CBG ( I )=C6G(.I )+SUM 
CALL DATS* (.3, JUNK) 

00 TU (• 3 2 0 ,' 3 3 U ) , JUNK 
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nnoo non 


320 WRITE ( I* P , 6 3 0 ) 

CALL AAA(.CBG,NP) 

330 I T R= ITR+ 1 

IF (,ITk- 30) ft 0 ,'8 0 ,*340 

NUN CONVERGENCE DUMP 

340 WRI TE< IP, 750) 

WRITE ( IP, 740) 

CALL AAA (XE LTA » NP ) 

WRH tl I-P ,730) 

CALL A A A (XV A L ,-NP ) 

DO 3b 0 J^Lt'NLVL 
NR-MAXX ( J ) 

!\ic=mxcol u ): 

T AUX= WHERE ( 2 , J ) 

WRITE! I P ,> 7-60 ) J , J , NR • NC , T AUX 
00 350 K=1,*NR 

350 w R I* T b ( IP, 7 70) (.ARK A Y C J , K , L ). , L = 1 , NC ). 

CALL EX I F 
360 DO 360 K = L ,-NP 

IF (. A L PH ( K ) - A L M AX (.K ) ) 360,370,370 

370 WRIXL ( I-P, 780 } K , ALPH ( K ) , ALMAX {.K ) 

IS-2 

IK r=i 

380 CONTINUE 

CALL DAI SW(.3, JUNK) 

GO )U (.390,400) , JUNK 
390 WRI- T L ( IP, >790) 

CALL A A A (. C L M A X , N P ). 

400 WRI I F( IP, 800) 

CALL AAA (.CVALt-NP ) 

DO A 10 K = 1,‘NP 

410 CL (. K ) = C L M A X ( K ).— C V A L ( K ) 

WRITE ( I P , o 1 0 ) 

CALL A A A (XL , NP ) 

SUN 1 = 0. 

DO 420 1^1,-NP 

420 SUM 1 = CBC (■ I )' * L T A ( I ). + SUM L 
0 = A S P E C * B W X * * 2 
COL I F r = U*SUM 1 

WRI* IE ( IP, -82 0) AL P HB,CULI FT 
GO TU (.10,430)., IS 
430 I AGN= 1 AGN+ 1 

IF (. IAGN-2) 4 50,4 5u, 440 
44 0 IF(l)lK) 460,4 50,4 70 

4 50 GO TU (.490,500), 1KT 
460 GO TO (.490,460), IKT 
470 GO TU (4 80, '900), IKT 
480 1 Q I = 2 

490 I Y= 1 

GO TO (500f95O),- iOV 
500 GO TU (510, -530, 530,530,550) , 1 A G N 

91.0 I P R 4 B = I P R A B - 1 

I F ( I PKAO-l ) 920,530,920 

520 A L PUB - AL PHB- ( (.ALPUV( 1 PR AL ) - A L PH V ( IPKAB-.l) )/2. ) 
GO TU 540 

5 30 A L P HO - (.ALP H 15 -» (,- 1 . )> * 1 K 1 *0 . 2 ) 

A L P H Li = A L P h B + D . 0 ( • 1 
I") I K = ( — 1 . • ) # * I K T 
540 wKl-rt ( i P, 960 ) N AML 
WRITE ( IP, 830) ALPh 0 
[ P R A b = I PRAB+ 1 
GO TU 2u 
5 50 RE I URN 


DID IT STALL . 
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noooonoooononononnnooono 


NO 


GE r 

ANOTHER 

ALPHA 


560 EUKMA T (. 1 h I / 1H0/35X , 25A2/ IX ) 

570 HGk^AT( 1 X * 40 ( 3H . / . )!/lH0/15X,30hbOUY ANGLL Of ATTACK, DEG. . =,F1C. 
12, LOX , 3i_HVALUF OF DISCRIMINANT. . . . = , E 10 . 6/ 1 5X , 30HB0D Y HL IGIHT / 

2SPAN. ..... =, F 10.2, IOX, 30HC0GY WIDTH / SPAN =,F10. 

32/15X,30HASPECI RATIO : . - , F 1 0 . 2 , 1 0 X , 3 OH W I N G HEIGHT / 

4SPAN = ,' F 1 0 . 2 / 1 5 X , 3 0 H w I N G BODY INCIDENCE, DEG . . =,F10 

5.2, IOX, 301 TIP THICKNESS CHORD. .... =,F1G.2) 

580 FORMAT ( 15X , 30HRU0T THICKNESS CHORD .... = , F 1 0 . 2 , 1 OX , 30HGEOKE TR I 
1C T vv 1ST, DEG .... = , F 10 . 2 / 1 5X , 30HNUMB t.l< OF SPaNWISE STATIONS. = 
2, FI 0.2, IOX, : 30H AERODYNAMIC TWIST, DEG . . . = , F 10 . 2/ 15X , 30HF LAP SPA 
3N / WING SPAN. . . . = , F 10 . 2 , 1 OX , 30HT APEK RATIO. = 

4 » F L 0 . 2 / 1 5 X ,. 3 0 H F L A P SETTING, DEG = , F 10 . 2 , 1UX , 30HRt YNOLDS 

5 NUMBER = * F 1 0 . 2 / 1 5 X *• 3 7 H C Gtl K D I N A T E S OF MOMENT REFERENC 

6F POINT, i.4X,2HX=,F10.2,9X,‘2HZ = ,F10.2) 

590 FORMAT ( 1 HO/ 1 X ,'40 !.3H . / . )7 IX 1 
600 FORMAT (. 1 OX ,-l 6HTH I CKNfcSS FACTOR) 

6 LO FORMAT MOX, 1 3 H R DO BAR / DU) 

620 FORMAT (.10X,-4HALPG ) 

630 FORMAT ( 1UX, -3HCBG) 

640 FOR HA T ( 1 OX ,'bHCLSTA ) 

650 FORMA V ( iOX , 2HA 3 ) 

660 FORMA T (.1 OX ,'5HCLSTU ). 

670 FORM AT ( 10X,'2HA4 ) 

680 FORMA I (. L0X. 2HA1 ) 

690 FORM ATT 10X.-2MA2) 

700 FORMA I ( 1GX,'4HALPH) 

710 F0RMAT<.1GX,'4HALPC) 

720 FORMAT MOX , 5HALPHU ) 

730 FORMAT MOX, 4HCVAL ) 

740 FORMAT! 1 OX ,'5HDt L T A ) 

750 FORMAT ( 1 X , 4LHUN ABLE TO CONVERGE AFTER 30 ITERATIONS ABORTED) 

760 FORMAT! 1HG/1H0,7H ARRAY ( , I 1 , 1 1 h 1, 1) TO ( , I 1 , 1 H , , I 2 , 1 H , , I 2 , 1H ) , F 

110.B/1H0) 

770 FORMAT! 121 10.4) 

760 FORMAT! IX, 1 IHSTALLEU AT ,I3,32H POINT SECTION ANGLE OF ATTACK =,F 
110.3, 29HM AX . SECTION ANGLE Of A I T ACK= , F 1 0 . 3 ) 

790 FUKMAT (. IOX , 5HC.LM AX ) 

600 FORMAT! IOX, '40HDISTRlbUT ION OF SECTION LIFT COEFFICIENT) 

810 FORMAT MOX, 25HS TALL MARGIN DISTRIBUTION) 

820 FORMAT MX, 120 MH. )/ lhO/ L5X, 30HBUDY ANGLE OF ATTACK, DEG. . =,F 10. 
12, 10X,30HLIFT COEFFICIENT ...... -»F10.6/1X,120MH. ) ) 

830 FORMA) !• I X , 4 0 ( 3 H . / .: ) / 1 X / 1 0 X ,■ 2 8 H A U J U S T E D ANGLE OF ATTACK TO ,F10.2, 

1 3 I HS LARCH I-NG FOR EXACT STALL PO I NT , / 1 X / 1 X , 4 0 ( 3H . / . ) / 1 X ) 

END 


YFS 

• 

BACK OFF 1/2 THE DISTANCE 
BETWEEN LAST ALPHA VALUE 
AND PRESENT ALPHA IF THIS 
WAS THE FIRST STALL 
OTHERWISE INCREASE OR 
DECREASE I HE ANGLE OF 
ATTACK BY 0.2 DEGREES 
DEPENDING ON. 

• 

IS IT STALLED AT I /Z POINT 

• • 

NC YtS 

SUBTRACT 0.2 ADD 0.2 

IF THIS IS THE FIRST ANGLE 
OF AT I AC K THEN DECREASE BY 
0.2 DEGREES UNTIL UNSTALL 
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APPENDIX B 


AVAILABILITY OF THE COMPUTER PROGRAM 


The computer program developed for the CDC 6400 and 6600 
series computers has been stored in the COSMIC Computer Center 
and is available for public use. All inquiries pertaining 
to availability and use of the program should be directed to: 

COSMIC 

Computer Center 
Barrow Hall 
University of Georgia 
Athens, Georgia 30601 

Telephone: (404) 542-3265 
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